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PREFACE 


Tne  CDC  computer  version  of  the  POEM  Monte  Carlo 
code  for  the  calculations  of  x-ray  photoemission  and  x-ray 
transition  zone  dose  is  now  available.  The  code  may  be 
obtained  by  sending  two  2400  foot  1/2  inch  reels  of  magnet- 
ic tape  to: 

General  Electric  Company  - TEMPO 
Attention:  DASIAC/ESPIG 

Post  Of  rice  Drawer  QQ 
Santa  Barbara 
California  93102 


The  tapes  must  be  accompanied  by  a letter  requesting  the 
code  and  a distribution  approval  letter  from  DNA  Headquar- 
ters. Approval  is  obtained  by  writing  to 

LTC  J.M.  Daley  (RAEV) 

Defense  Nuclear  Agency 

Washington,  D.C.  20305 

No  copy  of  the  POEM  code  will  be  released  until  DNA  approval 
has  been  obtained. 

The  development  of  the  POEM  code  is  continuing  - 
revisions  continue  to  be  incorporated.  Information  on  cod- 
ing revisions  will  be  distributed  through  SAI  Technical 
Notes.  Users  desiring  further  information,  other  versions 
of  the  code,  or  inclusion  in  the  Technical  Note  distribu- 
tion should  contact  the  author 

W.L.  Chadsey 

Science  Applications,  Inc. 

1651  Old  Meadow  Road 

McLean 

Virginia  22101 
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INTRODUCTION 


POEM  is  a FORTRAN  code  for  computing  x-ray  generated 
electron  currents  and  depositions  in  structures.  The  Monte 
Carlo  method  is  used  for  the  computation:  an  electron  source 

distribution  is  calculated;  a large  number  of  electron  tra- 
jectories are  then  traced  using  a random  walk  procedure,  and 
contributions  to  currents  and  depositions  are  scored  for 
each  trajectory.  After  a large  number  of  electron  histories, 
the  contributions  are  averaged  to  obtain  statistical  esti- 
mates and  standard  deviations  for  the  currents  and  depositions. 

The  primary  applications  of  the  POEM  code  have  been 
the  calculation  of  x-ray  photoemission  from  surfaces  for 
IEMP  (internal  electromagnetic  pulse)  and  SGEMP  (system- 
generated electromagnetic  pulse)  analyses  and  the  calcula- 
tion of  transition  zone  dose  and  current  for  TREE  (transient 
radiation  effects  on  electronics)  analyses.  These  are  the 
applications  for  which  the  code  was  originally  developed 
and  for  which  most  interest  has  been  expressed.  The  pro- 
cedures for  these  applications  of  the  code  are  discussed  in 
detail  in  this  report.  Other  applications  of  the  POEM  code 
have  included  analyses  of  x-ray,  y-ray,  and  SEM  (scanning 
electron  microscope)  produced  electron  currents  and  deposi- 
tions in  multilayer  devices,  transient  currents  and  fields 
in  x-ray  and  electron  irradiated  dielectrics,  induced  cur- 
rents on  conductors  in  x-ray  and  y-ray  irradiated  plated 
wire  memories,  and  x-ray  and  y-ray  produced  currents  and 
depositions  in  air  near  conductor  surfaces.  These  applications 


8 


are  only  briefly  discussed  in  this  report;  additional  infor- 
mation can  be  obtained  from  the  cited  references  and  from 
the  author. 

The  POEM  code,  as  stated  in  the  title  of  this 
report,  is  a "fast"  Monte  Carlo  code.  This  requires  some 
explanation.  Over  the  past  several  years,  a number  of  trans- 
port codes  have  been  developed  which  are  applicable,  for 
example,  to  the  calculation  of  x-ray  photoemission.  These 
include  Monte  Carlo  codes  - POEM,'*'  SANDYL,^  ETRAN,^  and 
FASTER-3ETA4  - and  analytical  codes  - QUICKE2,5  EASY1,6 
and  BOLT. ' The  Monte  Carlo  codes  are  the  most  rigorous, 
accurate,  and  widely  explicable.  Yet,  Monte  Carlo  codes 
are  in  general  large,  complicated  codes  which  require  large 
amounts  of  computer  core  and  time  to  produce  predictions 
containing  statistical  uncertaintly . Recently,  research 
has  concentrated  on  developing  analytical  codes  which  offer 
the  appeal  of  the  more  elegant  solution,  reduced  computation 
cost,  and  improved  precision.  Of  necessity,  however,  the 
analytical  codes  rely  on  a number  of  simplifying  approxima- 
tions. At  best,  these  approximations  restrict  the  applica- 
bility of  the  codes;  at  worst,  they  invalidate  the  predic- 
tions through  oversimplification.  Early  codes  such  as 
8 9 

GRAP  and  EGRESS  rely  on  an  invalid  rectilinear  electron 
transport  model  which 'produces  erroneous  predictions  for 
back  emission  yield  and  fore  and  back  emission  angular  dis- 
tributions. While  these  codes  have  been  known  to  be  in 

10 

error  for  some  time,  their  use  unfortunately  remains 
common.  Later,  more  sophisticated  codes  use  a truncated 
Legendre  series  to  represent  the  electron  distribution 
function.  The  EASY1  code  of  Lopez  uses  a Pi  approximation; 
the  QUICKE2  code  of  Dellin  and  MacCallum  uses  a P2  approxi- 
mation. For  the  conditions  for  which  the  code  was  designed 
QUICKE2  is  an  accurate  and  exceedingly  convenient  tool. 
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It  must  be  kept  in  mind  thouqh  that  its  applicability  is 
limited.  QUICKE2  predictions  disagree  with  Monte  Carlo 
predictions  and  experiment  for  back  photoemission  yield 
and  fore  and  back  emission  angular  distributions  for  photon 
energies  greater  than  100  keV.^  More  important,  none  of  the 
existing  analytical  codes  provide  information  on  the  two- 
dimensional  angular  distribution  of  emission  for  the  case 
of  non-normal  photon  incidence,  and  none  provide  information 
on  the  dose  distribution  profile  near  a material  interface. 

Rather  than  abandon  the  Monte  Carlo  method  and 

completely  embrace  the  convenient,  but  limited  analytical 

codes,  an  alternative  course  has  been  pursued.  We  have 

built  a fast  Monte  Carlo  code  for  the  calculation  of  x-ray 

photoemission  and  transition  zone  dose  and  current.  We 

have  optimized  the  POEM  Monte  Carlo  code  and  incorporated 

variance  reduction  techniques  such  as  stratified  sampling, 

importance  sampling,  track  splitting,  and  control  variates 

which  have  resulted  in  an  order  of  magnitude  reduction  in 

computation  time  requirements.  For  example,  for  the  POEM 

code  calculations  of  x-ray  photoemission  presented  by 

2 2 

Garth  and  Chadsey,-  m which  comparisons  were  made  with 

13 

tne  measurements  of  Bradford,  approximately  15  CP  minutes 
were  required  on  a CDC  6 600  to  obtain  a one  percent  to  two 
percent  standard  deviation  on  the  prediction  of  total  yield 
and  five  percent  to  ten  percent  on  the  prediction  of  the 
emission  distributions.  This  same  precision  can  now  be 
obtained  with  approximately  one  CP  minute  on  the  CDC  6600. 

This  dramatic  increase  in  code  efficiency  has  al- 
lowed us  to  perform  extensive  parametric  studies  of  x-ray 
photoemission  and  transition  zone  dose  and  current.  The 
emission  study  included  calculations  of  fore  and  back 
emission  yield,  energy  distribution,  one-  and  two-dimen- 
sional angular  distributions,  and  the  various  joint 
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energy-angular  distributions  for  the  elements  C,  A1 , Ti,  Cu, 
Mo,  Sn,  Ta,  and  Au,  photon  energies  0.01,  0.02,  0.05,  0.10, 
0.20,  0.50,  1.0,  2.0  MeV  and  E *5  keV  where  Ev  is  the  K-edge, 

I\  I\ 

and  photon  angles  of  incidence  0,  30,  60,  and  80  degrees 
with  respect  to  the  surface  normal.  The  study  also  included 
comparisons  of  the  POEM  calculations  with  QUICKE2  calculations 
over  the  range  of  the  parametric  study. 

The  transition  zone  dose  and  current  parametric 
study  was  performed  for  gold/polyethylene  and  gold/silicon 
interfaces  for  irradiation  by  monochromatic  photon  spectra 
in  the  range  0.01  to  2.0  Me\7.  The  results  include  dose  and 
current  profiles,  contributions  due  to  electron  emission 
and  backscatter  at  the  interface,  contributions  due  to  the 
various  source  electron  groups  - K and  L photoelectrons,  K 
Auger  electrons,  and  Compton  electrons  - and  analytical 
functional  fits  to  the  dose  and  current  profiles. 

Extensive  computational  studies  such  as  the  photo- 
emission and  transition  2one  dose  and  current  parametric 
studies  would  not  have  been  feasible  with  the  old  version 
of  the  POEM  code  or  any  of  the  other  transport  codes,  ETRAN, 
SANDYL,  or  FASTER-BETA.  Computation  time  requirements 
would  have  been  excessive.  Thus,  the  conversion  of  the  POEM 
code  to  a fast  Monte  Carlo  code  has  greatly  increased  its 
power  as  a computational  tool.  To  be  fair,  however,  it  must 
be  pointed  out  that  the  increase  in  code  efficiency  could 
only  be  obtained  through  some  sacrifice  in  generality. 

Codes  like  SANDYL  and  FASTER-BETA  are  large  codes  applicable 
to  general  types  of  problems.  If  one  wants  to  make 
frequent  calculations  for  a problem  for  which  the  POEM 
code  was  designed,  one  would  be  advised  to  use  POEM  because 
of  the  significant  computation  cost  savings.  On  the  other 
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hand,  for  infrequent  calculations  for  complex  geometry 
problems  otner  than  those  for  which  POEM  was  designed,  one 
would  be  advised  to  use  a code  like  SANDYL.  While  we  have 
had  little  difficulty  in  adapting  POEM  to  various  complex 
geometries,  such  as  the  plated  wire  memory,  this  requires 
familiarity  with  the  code  and  of  course  an  investment  of 
labor . 

Three  reports  on  the  POEM  code  are  being  published. 
This  report  presents  the  code  documentation.  The  other 
two  reports  present  the  results  of  the  x-ray  photoemission 
and  transition  zone  dose  and  current  parametric  studies. 

The  documentation  and  release  of  the  POEM  code  should  pro- 
vide the  technical  community  with  a convenient  and  accurate 
tool  for  IEMP,  SGEMP,  and  TREE  analyses.  The  results  of 
the  x-ray  photoemission  study  should  serve  to  acquaint 
researchers  with  the  characteristics  of  electron  emission 
back  and  fore  yields,  energy  spectra,  and  one-  and  two- 
dimensional  angular  distributions  as  functions  of  material 
atomic  number  and  photon  energy  and  angle  of  incidence,  and 
to  indicate  the  ranges  of  applicability  of  the  analytical 
codes.  The  results  of  the  transition  zone  dose  and  current 
study  should  serve  to  acquaint  researchers  with  the  charac- 
teristics of  dose  and  current  profiles  near  high-Z  low-Z 
interfaces  and  to  serve  as  a data  base  against  which  to 
test  future  analytical  or  empirical  models. 

The  POEM  code  started  several  years  ago  as  a version 
of  ETRAN  which  we  specialized  for  the  calculation  of  x-ray 
photoerr.ission  and  transition  zone  dose  and  current  by  delet- 
ing all  procedures  extraneous  to  these  calculations.  The  code 
has  undergone  continuous  revision  since  that  time.  The  elec- 
tron source  routines  were  replaced,  the  electron  multiple 
scattering  distribution  calculation  was  revised,  and  the 
random  walk,  scoring,  and  output  routines  were  completely 
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rewritten.  Recently,  the  code  was  optimized,  requiring 
rewrite  of  essentially  the  entire  code,  variance  reduction 
techniques  were  incorporated,  and  fitting  routines  were 
added  to  the  transition  zone  dose  and  current  calculations. 
Revision  of  the  code  has  by  no  means  been  completed;  in 
fact,  at  this  time  we  are  working  on  a revolutionary  change 
in  the  Monte  Carlo  approach.  Documentation  and  release  of 
the  POEM  code,  however,  is  long  overdue.  This  report, 
therefore,  describes  POEM  as  of  January  1975  - the  version 
used  in  the  parametric  studies.  In  the  future,  modifica- 
tions to  the  code  will  be  described  in  Technical  Notes  to 
be  distributed  to  the  users. 

Section  2 of  this  report  describes  the  code  appli- 
cations - problem  geometry,  range  of  validity,  and  computed 
quantities.  Section  3 describes  the  computational  proce- 
dures - source  generation,  history  origination,  random  walk, 
scoring,  history  termination,  statistical  estimation  of 
output  quantities  and  standard  deviations,  and  variance 
reduction.  Section  4 describes  the  photon  and  electron 
cross  sections  used.  Section  5 contains  the  user  instruc- 
tions and  input/output  description.  The  cited  references 
may  be  found  at  the  end  of  Section  5.  Code  listings  and 
sample  problems  are  contained  on  the  master  tape  along 
with  the  program  source  deck  images. 

POEM  is  a family  of  electron  transport  codes  con- 
tained in  an  UPDATE  program  library.  Each  code  in  the 
library  is  designed  for  a specific  application:  x-ray 

photoemission,  x-ray  transition  zone  dose,  SEM  dose,  dose 
and  charge  deposition  in  dielectrics,  radiation-induced 
currents  in  plated  wire  memories,  and  so  on. 

The  codes  are  filed  by  UPDATE  deck  name.  The  deck 
name  for  each  application  is  specified  in  Section  2. 
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Complete  code  descriptions  are  included  in  Section  2 for 
the  x-ray  photoemission  and  x-ray  transition  zone  dose  ap- 
plications. Brief  code  descriptions  are  also  included  in 
Section  2 for  a number  of  other  applications.  Additional 
information  on  these  codes  may  be  obtained  by  contacting 
the  author. 
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Section  2 


APPLICATIONS 


2.1  X-RAY  PHOTOEMISSION:  (DECK  NAMES  1S1D,  1S2D) 

X-rays  and  y-rays  produce  swift  Compton,  photo-,  and 
Auger  electrons  in  matter.  These  primary  electrons  in  turn 
produce  secondary,  knock-on  electrons.  Those  electrons  pro- 
duced within  an  electron  range  from  a material  surface  have 
a non-zero  probability  of  escaping  from  that  surface.  These 
x-ray  and  y~ray  generated  electron  emissions  drive  IEMP, 

SGEMP , and  source  region  EMP.  Accurate  prediction  of  the 
IEMP,  SGEMP,  or  source  region  EMP  fields  depends  on  the  ac- 
curate determination  of  the  yield  and  energy-angular  distri- 
bution of  the  electron  emissions.  Prediction  of  electron 
emission  requires  calculation  of  the  electron  source  in  the 
region  within  an  electron  range  from  the  emitting  surface 
and  the  calculation  of  the  electron  transport  to  the  surface. 
POEM  uses  an  analytical  model  for  the  source  calculation 
and  a Monte  Carlo  model  for  the  transport  calculation. 


2 1.3 


Problem  Geometry 


POEM  calculates  the  electron  emission  from  a planar 
surface.  The  restriction  to  planar  geometry  is  not  a severe 
limitation.  If  the  radius  of  curvature  of  a surface  is 
large  compared  with  the  electron  range,  the  surface  effec- 
tively appears  planar  in  the  transport  calculation.  Thus, 
the  code  can  be  used  to  calculate  the  electron  emission  for 
such  geometries  as  spheres,  cylinders,  and  conical  sections 
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provided  the  radius  of  curvature  is  large  compared  with 
an  electron  range.  Since,  for  example,  electron  ranges  in 
aluminum  for  0.01,  0.10,  and  1.0  MeV  electrons  are  0.00013, 
0.0069,  and  0.20  cm,  respectively,  this  is  not  a particularly 
restrictive  condition.  One  exception  is  the  electron  emis- 
sion from  a small  diameter  wire.  For  this  case,  a specie] 


version  of  the  code  has  been  built  as  described  in  Section  2.3 


The  problem  geometry  for  fore  and  back  emission  is 
defined  in  Figure  1.  The  emitting  surface  of  a planar 
slab  of  irradiated  material  lies  in  the  xy-plane;  the  z-axis 
lies  on  the  outward  normal  to  the  surface.  The  photons  are 


incident  parallel  to  the  xz-plane  with  angle,  0q,  with  re- 


spect to  the  z-axis  and  with  positive  direction  cosine 
with  respect  to  the  x-axis.  POEM  calculates  the  electron 
emission  for  the  case  of  plane  wave  photon  irradiation. 
(For  narrow  beam  photon  irradiation,  the  POEM  calculated 


yields  should  be  multiplied  by  sec0Q.)  The  photon  flux  is 


assumed  to  be  defined  at  the  emitting  surface.  The  POEM 
code  does  not  calculate  the  photon  flux  attenuation  within 
the  irradiated  material.  (The  flux  attenuation  can  be 


easily  calculated  using  any  of  a number  of  readily  avail- 

31  32 

able  photon  transport  codes  such  as  FSCATT  or  MORSE.  ) 


Since  the  photon  flux  is  assumed  to  be  defined  at  the  emit- 
ting surface,  the  electron  emission  yield  is  independent  of 
the  slab  thickness  for  thicknesses  greater  than  the  maximum 
electron  range.  Unless  specified,  the  slab  thickness  is 
assumed  equal  to  the  maximum  electron  range.  For  slab 
thicknesses  thinner  than  the  maximum  electron  range,  the 
thickness  must  be  specified. 

An  electron  is  emitted  in  direction  & with  angle  0 
with  respect  to  the  z-axis  and  angle  ({>  with  respect  to  the 


xz-plane.  Emission  from  the  slab  surface  nearer  the  source 
is  termed  ba^k  emission  (Figure  lb) . Emission  from  the 
slab  surface  farther  from  the  source  is  termed  fore  emission 
(Figure  la)  . 

It  is  important  that  the  POEM  code  calculates  the 
full  two-dimensional  angular  distribution  of  yield  y'(0,4>).* 
The  analytical  codes  EASil  and  QUICKE2  calculate  only  the 
one-dimensional  angular  distribution  y'(0)  which  is  the  full 
angular  distribution  integrated  over  azimuth  <J>.  The  analy- 
tical codes  provide  no  information  on  the  azimuthal  dis- 
tribution of  emission  for  non-normal  photon  incidence. 
Knowledge  of  the  full  two-dimensional  angular  distribution 
is  critical  to  practical  IEMP,  SGEMP,  and  source  region  EMP 
problems.  Consider  for  example,  the  SGEMP  for  a right- 
circular  cylindrical  structure  irradiated  normal  to  the 
cylindrical  axis  by  plane  wave  x radiation  (Figure  2) . 

At  any  point  of  the  irradiated  cylindrical  surface,  except 
along  the  line  tangential  to  the  wave  front,  the  electron 
emission  angular  distribution  is  asymmetric  about  the  sur- 
face normal  at  that  point.  At  each  such  point,  there  is 
therefore  a component  of  the  emission  current  tangential  to 
tr.e  cylindrical  surface  as  well  as  normal  to  the  surface. 
This  component  drives  circumferential  surface  currents  on 
the  structure.  The  analytical  codes  cannot  predict  the 
magnitude  of  this  tangential  component  of  the  emission 
current  since  they  provide  no  information  on  the  azimuthaJ 
angular  distribution  of  electron  emission.  The  POEM  code', 
on  the  other  hand,  does  predict  both  the  normal  and  tan- 
gential components  of  the  emission  current. 


* Notation:  y denotes  the  total  electron  yield  (electrons/ 
photon);  y'  denotes  a distribution  of  the  yield  with  respect 
to  one  or  more  variables. 
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2.1.2 


Limitations 


Photon  Environment  - The  incident  photon  spectrum 
can  be  specified  in  one  of  three  ways:  (1)  monochromatic 

spectrum,  (2)  black  body  spectrum,  or  (3)  an  arbitrary 
spectrum  defined  by  up  to  120  photon  energy  groups.  The 
high  energy  limit  for  the  photon  spectrum  is  10  MeV.  This 
limitation  is  principally  imposed  by  the  fit  used  for  the 
radiative  energy  loss.  POEM  does  not  include  pair  produc- 
tion in  the  electron  source  calculation.  Because  of  the 
positron/electron  symmetry  of  production  and  near  symmetry 
of  transport,  the  pair  production  contribution  to  the  emis- 
sion current  is  insignificant.  For  predictions  for  an 
emission  experiment  in  which  positrons  and  electrons  are 
distinguishable,  however,  care  must  be  exercised  in  using 
POEM  for  photon  energies  greater  than  several  MeV,  for  which 
energies  pair  production  constitutes  a sizeable  fraction  of 
the  total  photon  interaction  cross  section.  The  low  energy 
limit  for  the  photon  spectrum  is  determined  by  the  low  energy 
limit  for  the  electron  transport.  The  low  energy  cutoff  on 
the  electron  transport  is  currently  1 keV.  A practical  low 
energy  limit  on  the  photon  spectrum  is,  then,  about  5 keV. 
While  the  electron  scattering  cross  sections  used  are  not 
strictly  valid  at  electron  energies  as  low  as  1 keV,  in 
practice  we  have  found  reasonable  agreement  of  prediction 
with  experiment  for  photon  energies  on  the  order  of  10  keV. 


As  discussed  previously,  the  photon  spectrum  is 
assumed  to  be  defined  at  the  electron  emitting  surface. 
For  cases  in  which  there  is  a significant  thickness  of 
material  between  the  photon  source  and  the  electron  emit- 
ting surface,  the  photon  intensity  and  spectrum  at  the 
surface  must  be  computed  using  an  appropriate  photon 
transport  code. 
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POEM  assumes  plane  wave  photon  radiation  incident 
on  the  electron  emitting  surface.  The  POEM  predictions, 
however,  are  also  applicable  to  narrow  beam  photon  irradi- 
ation of  the  surface  provided  that  the  calculated  electron 
yields  are  multiplied  by  sec0Q,  where  6^  is  the  photon 
angle  of  incidence  with  respect  to  the  surface  normal. 

The  photon  angle  of  incidence  with  respect  to  the 
surface  normal  is  arbitrary,  but  the  incidence  must  be 
monodirectional . Extension  of  the  code  to  accept  an  angu- 
lar distribution  of  photon  incidence  requires  only  minor 
coding  changes  provided  that  the  photon  angular  distribu- 
tion and  energy  distribution  can  be  assumed  to  be  uncor- 
related. Extension  to  the  case  of  correlated  photon  angu- 
lar and  energy  distributions  would  require  major  code 
modification. 

POEM  computes  both  fore  and  back  electron  emission 
for  a given  photon  environment,  or  if  so  requested  by  the 
user,  only  fore  or  back  emission. 

Material  Configuration  - The  electron  emitting 
material  may  be  any  homogeneous  material  composed  of  a 
single  element  or  a composition  of  up  to  ten  elements. 

The  elements  may  be  any  of  those  for  which  the  photoelectric 
cross  sections  are  compiled  on  the  ENDF/B  file:  these  are 

the  87  elements  with  atomic  number  Z=1  through  83,  86,  90, 
92,  and  94. 

POEM  assumes  the  electron  emitting  surface  to  be 
planar.  The  surface  geometry  may  be  multidimensional,  as 
discussed  previously,  provided  that  the  radius  of  curvature 
of  the  surface  is  large  everywhere  compared  with  the  maxi- 
mum electron  range. 
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The  electron  emission  may  be  calculated  for  a homo- 
geneous slab  of  material  of  arbitrary  thickness.  If  the 
slab  is  thick  compared  with  the  maximum  electron  range, 
the  electron  yield  is  independent  of  thickness,  since  the 
photon  intensity  and  spectrum  are  defined  at  the  emitting 
surface.  The  emission  may  also  be  calculated  for  a foil 
tninner  than  the  maximum  electron  range,  provided  that  the 
transport  calculation  is  set  up  so  that  the  maximum  step- 
size  in  the  random  walk  is  reasonably  small  compared  with 
the  material  thickness.  If  the  emitting  structure  is  con- 
structed of  several  layers  of  different  materials,  the 
surface  layer  of  which  is  thinner  than  the  electron  range, 
then  the  emission  must  be  calculated  using  the  multislab 
version  of  the  code  (Section  2.3). 


2.1.3  Computed  Quantities 

The  x-ray  photoemission  quantities  computed  by  POEM 
are  summarized  in  Table  1.  All  quantities  shown  in  this 
table  are  computed  for  both  fore  and  back  emission.  The 
total  yield  y is  expressed  in  electrons/photon,  coul/cal (pho- 
tons) and  cal  (electrons)  /cal  (photons)  . If  y'(E,0,<|>)  is  the 
distribution  of  yield  per  unit  electron  energy,  per  unit 
solid  angle,  per  incident  photon,  the  calculated  total 
yield  is 


Jraax  2tt  tt/2 


y = y dE J dtyj'  dO  sinO  y ' (E  , 0 , 


*) 


nun 


(1) 


where  Em^n  = 1 keV  is  the  minimum  electron  energy  and 

E^„  ilhv  is  the  maximum  electron  energy  and  hv  is  the 
max  max  1 max 

maximum  photon  energy. 


22 


Table  1 (Continued) 


Table  1 (Continued) 
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probability  of 
from  source 


The  various  emission  distribution  functions  are  ex- 
pressed both  in  absolute  units  and  normalized  to  the  total 
yield.  The  one-dimensional  distribution  functions  are 
defined  as 

2tt  tt/2 


y'  (e) 


■/d*/ 


d0  sin6  y ' (E,  0 , <j>) 


(2) 


y ' (0) 


2tt 


E 


max 


dE  y ' (E,  6 , <{>) 


min 


(3) 


y'  (<l>) 


tt/2  1 

- f de/ 


max 

dE  y' (E , 6 , <j>) 


mm 


(4) 


Note  that  the  polar  angular  distribution  is  expressed  both 
as  the  electron  yi_id  per  unit  solid  angle  y1  (0)  and  as  the 
electron  yield  per  unit  scattering  angle  y"(0)  = 2tt  sinO  y'(0 

The  two-dimensional  distribution  functions  are  de- 
fined as 


2tt 

y’ (E,6)  = f d<£>  y 1 { E , 0 , <J> ) 
o 


y'  (E,<p) 


sin0  y ' (E , O , <|>) 


o 


(5) 


(0) 
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The  normal  maximum  resolution  for  the  electron 
emission  distribution  functions  is  20  energy  bins,  10  9- 
angle  bins,  and  8 <J>-angle  bins.  Resolution  can  be  increas- 
ed through  increasing  the  dimensions  of  the  code,  as  for  ex- 
ample was  done  in  the  calculations  reported  by  Garth  and 

12 

Chadsey . 

Several  moments  of  the  emission  distribution  are 
also  calculated.  These  include  the  mean  emitted  electron 
energy  <E>,  the  mean  cosines  of  the  polar  angular  and  azi- 
muthal angular  distributions,  <cos0>  and  <co:;<j>>,  and 
<sin0cos(f)> . The  <cos9>  is  a measure  of  ho*,  strongly  the 
emission  solid  angular  distribution  is  directed  in  the 
direction  of  the  outward  normal  to  the  surface.  The 
<cos<f>>  and  <sinbcos<f>>  are  measures  of  the  asymmetry  of  the 
azimuthal  angular  distribution  due  to  oblique  photon  inci- 
dence. (For  normal  photon  incidence  <cosc{>>  and  <sin0cos<})> 
are  zero . ) 

Since  POEM  is  a Monte  Carle  code,  the  computed 
emission  quantities  are  statistical  estimates  of  the  emis- 
sion parameters.  It  is  therefore  important  not  only  to 
calculate  the  desired  quantities  bur.  also  to  estimate  tne 
statistical  precision  of  the  calculations.  For  each  of 
the  electron  emission  quantities  discussed  above,  POEM 
calculates  the  estimated  standard  deviation. 

In  addition  to  these  electron  emission  quantities, 
POEM  calculates  several  quantities  describing  the  electron 
source  and  the  electron  transport.  The  source  parameters 
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are  the  source  density  S in  electrons/photon  (averaged  over 
the  source  region),  the  source  energy  spectrum  S'(E),  the 
source  angular  distribution  with  respect  to  the  photon 
direction  S'(0),  and  the  two-dimensional  source  distribution 
S'(E,0).  Also  calculated  are  the  mean  source  electron 
energy  and  the  maximum  energies  of  the  Compton  electrons, 

K and  L photoelectrons,  and  K-Auger  electrons.  The  trans- 
port parameters  are  the  average  number  of  random  walk  steps 
and  average  total  deflection  of  the  particles  which  escape 
the  surface,  the  total  number  of  particles  which  escape 
the  surface,  and  the  mean  probability  of  escape  for  the 
source  electrons. 

Two  versions  of  the  emission  code  exist:  one  for 

the  calculation  of  the  two-dimensional  angular  distribution 
y'(6/^)  (deck  name  1S2D) ; the  other  for  only  a one-dimensional 
calculation  y'(6)  (deck  name  1S1D) . For  cases  of  normal 
photon  incidence  or  for  which  the  two-dimensionality  of  the 
i emission  distribution  is  not  of  concern,  the  one-dimensional 

l 

version  should  be  used  since  it  is  significantly  faster. 

2.2  TRANSITION  ZONE  DOSE  AND  CURRENT:  (DECK  NAME  2S) 

Dose  is  defined  by  the  International  Commission  on 
Radiological  Units  (ICRU)  as  the  energy  imparted  to  matter 
by  ionizing  radiation  per  unit  mass  of  matter.  The  "energy 
imparted  to  matter"  is  that  which  appears  as  ionization  or 
i excitation,  increase  in  crystal-lattice  energy,  and  so 

| forth  in  the  material.  This  is  just  that  energy  removed 

I from  the  radiation  field  as  a whole.  For  y-ray  or  x-ray 

i irradiation  of  matter,  the  radiation  field  consists  of  the 

primary  photons  and  secondary  photons,  and  the  photo-Compton 
; and  secondary  electrons.  In  a region  of  a homogeneous 

material  far  from  the  photon  source  and  farther  from  any 
interface  with  a dissimilar  material  than  the  range  of  the 
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most  penetrating  of  the  electrons,  electron  equilibrium 
exists:  the  energy  transported  into  an  incremental  volume 

of  matter  by  electrons  is  equal  m the  mean  to  the  energy 
transported  out  of  the  volume  by  electrons.  The  equilib- 
rium dose  is  thus  equal  to  the  energy  removed  from  the 
photon  radiation  field  per  unit  mass  of  matter  - the  kerma. 
The  equilibrium  dose  can  be  determined  directly  from  cal- 
culation of  the  photon  transport;  the  electron  transport 
need  not  be  calculated. 

Near  a boundary  between  dissimilar  materials  elect- 
ron equilibrium  fails  because  of  differences  in  the  electron 
production  and  transport  properties  of  the  different  mate- 
rials. The  energy  transported  into  an  incremental  volume 
of  matter  by  electrons  is  not  equal  to  the  energy  trans- 
ported out  by  electrons.  Thus  the  dose  is  not  equal  to 
the  kerma.  Since  in  general  the  energy  directly  imparted 
to  the  matter  by  the  primary  photons  is  negligible  compared 
with  the  energy  imparted  by  the  electrons,  the  dose  is 
critically  sensitive  to  the  electron  balance.  Nonequilib- 
rium dose  can  differ  from  the  kerma  by  more  than  an  order 
of  magnitude. 

In  a region  of  a x-ray  or  y-xay  irradiated  material 
in  which  electron  equilibrium  exists,  the  photo-Compton 
current  is  proportional  to  the  photon  flux.  The  constant 
of  proportionality  is  a function  of  the  photon  energy  and 
the  atomic  number  of  the  material,  but  is  independent  of 
the  material  density  or  spatial  location  within  the  material. 
Near  the  boundary  between  dissimilar  materials,  this  rela- 
tionship breaks  down  because  of  the  failure  of  electron 
equilibrium.  The  nonequilibrium  photo-Compton  current, 
just  as  the  dose,  may  differ  significantly  from  the  equilib- 
rium value.  The  divergence  in  electron  current  near  a 
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boundary  between  dissimilar  materials  results  in  charge 
deposition  near  the  boundary. 

The  region  near  a material  interface  in  which  elec- 
tron equilibrium  fails  is  known  as  the  transition  zone. 

The  dose  and  current  distributions  in  the  transition  zone 
are  referred  to  variously  as  transition  zone  dose  and  cur- 
rent, non-local  charge  and  energy  deposition,  and  dose 
enhancement.  Transition  zone  dose  and  current  effects  are 
important  in  weapons  radiation  effects  problems  such  as 
TREE,  IEMP  in  cables,  and  charge  buildup  in  dielectrics, 
as  well  as  in  problems  in  radiography  and  radiology.  The 
POEM  code  provides  a convenient,  accurate  tool  for  the 
analysis  of  transition  zone  dose  and  current. 

2.2.1  Problem  Geometry 

POEM  calculates  the  transition  zone  dose  and  cur- 
rent in  the  region  within  an  electron  range  from  planar 
interface  between  two  materials  irradiated  by  plane  wave 
x or  y radiation.  The  problem  geometry  is  defined  in 
Figure  3.  The  planar  interface  between  materials  A and 
B lies  in  the  xy-plane.  The  z-axis  is  directed  away  from 
material  A.  The  plane  wave  photon  radiation  is  incident 
at  an  angle  0Q  with  respect  to  the  z-axis.  The  photon 
incidence  can  be  either  through  material  A or  B.  Just  as 
with  the  emission  code,  the  photon  intensity  and  spectrum 
are  assumed  to  be  defined  at  the  interface.  A and  B are 
planar  slabs  of  arbitrary  thickness.  N deposition  zones  of 
of  thickness  Az  are  defined  in  B by  planes  parallel  to  the 
interface  at  distances  Z\  = Az,  Z2  = 2Az,  . ..,  ZN  = NAz. 

The  photo-Compton  charge  and  energy  currents  are  calculated 
at  the  planes  z=0,  z\,  z2,  ...,  z^.  The  charge  and 
energy  deposition  are  calculated  in  the  N zones. 
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2.2.2 


Limi tations 


Photon  Environment  - The  transition  zone  dose  and 
current  calculation  is  valid  for  photon  energies  in  the 
range  5 keV<hv< 2 MeV.  The  lower  limit  is  determined  by 
the  1 keV  cutoff  on  the  electron  transport.  The  upper 
limit  is  determined  by  the  exclusion  of  pair  production  in 
the  electron  source  calculation.  While  pair  production 
does  not  significantly  contribute  to  the  photo-Compton 
current  because  of  particle/antiparticle  symmetry,  pair 
production  does  significantly  contribute  to  dose  for 
h^-2  MeV. 

The  calculation  assumes  plane  wave  photon  radiation. 
The  angle  of  incidence  with  respect  to  the  normal  to  the 
interface  can  be  CKO  <1°0<>.  Just  as  in  the  emission  calcu- 
lation, the  code  can  easily  be  modified  to  accept  an  angular 
distribution  of  photon  incidence,  provided  that  the  photon 
energy  and  angular  distributions  can  be  assumed  to  be  un- 
correlated. 

The  photon  intensity  and  spectrum  are  assumed  to 
be  defined  at  the  material  interface.  If  these  are  un- 
known, then  the  flux  attenuation  between  the  photon  source 
and  the  interface  must  be  calculated  using  an  appropriate 
photon  transport  code. 

Material  Configuration  - The  code  calculates  the 
transition  zone  dose  and  current  in  the  neighborhood  of  a 
planar  interface  between  two  materials.  The  materials  may 
be  any  of  those  materials  for  which  the  emission  code  can 
be  used.  The  code  normally  sets  the  material  thicknesses 
equal  to  the  maximum  electron  range  in  each  material. 

Since  the  photon  flux  is  assumed  to  be  defined  at  the 
interface,  increasing  the  material  thicknesses  beyond  the 
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maximum  electron  ranges  does  not  affect  the  results.  The 
minimum  material  thicknesses  are  determined  by  the  step 
sizes  selected  for  the  transport  calculation.  The  deposi- 
tion zone  width  should  be  larger  than  the  longest  step  size 
of  the  random  walk. 

2.2.3  Computed  Quantities 

The  calculated  transition  zone  dose  and  current 
quantities  are  summarized  in  Table  2.  The  current  and 
deposition  distributions  are  calculated  in  material  B as 
functions  of  the  distance  z from  material  A,  over  the 
region  from  the  interface  out  to  the  maximum  electron 
range.  The  distributions  are  the  photo-Compton  current 
per  unit  photon  flux  (electrons/photon) , the  photo-Compton 
energy  flux  per  unit  photon  energy  flux  (dimensionless) , 
the  charge  deposition  density  per  unit  photon  flux  (elect- 
rons/g) / (photon/cm2) , the  dose  per  unit  photon  energy  fluence 
(cm2/g) , and  the  relative  dose  (tbo  ratio  of  dose  to  kerma , 
which  is  dimensionless) . 

In  addition  to  the  dose  and  current  profiles,  the 
electron  yield  y at  the  interface  and  the  energy  distribu- 
tion y'(E)  and  angular  distribution  y'(0)  with  respect  to 
the  interface  normal  are  calculated.  For  the  case  of  the 
electrons  arising  in  the  deposition  material  (material  B) , 
the  yield  y is  the  backscatter  yield.  For  the  case  of  the 
electrons  arising  in  the  neighboring  material  (material  A) , 
the  yield  y is  the  emission  yield. 

Separate  calculations  are  made  of  the  contributions 
to  the  dose  and  current  profiles  due  to  (1)  emission  of 
electrons  across  the  interface  and  (2)  backscatter  of 
electrons  from  the  interface.  Separation  of  the  emission 
and  backscatter  components  of  the  profiles  is  useful  to 
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Table  2 COMPUTED  QUANTITIES:  X-RAY  TRANSITION  ZONE  DOSE 
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Table  2 (Continued) 
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understanding  the  interface  effects.  Furthermore,  the 
efficiency  of  the  Monte  Carlo  calculation  is  increased 
through  the  separation  of  the  components  and  the  proper 
allocation  of  the  Monte  Carlo  histories  to  the  separate 
calculations,  as  discussed  in  Section  3.2. 

Separate  calculations  can  also  be  made  of  the  con- 
tributions to  the  current  and  dose  profiles  due  to  electrons 
from  each  of  the  source  groups:  Compton,  K-photo,  L-photo, 

and  K-Auger  electrons.  This  separation  is  again  useful  to 
the  interpretation  of  the  interface  effects.  It  can  also 
be  used  to  increase  code  efficiency  and  to  increase  the 
resolution  of  the  dose  and  current  profiles. 

The  Monte  Carlo  output  is  in  the  form  of  the  charge 
and  energy  currents  at  each  plane  bounding  the  deposition 
zones  and  charge  and  energy  deposition  within  each  deposi- 
tion zone.  In  addition  to  the  Monte  Carlo  data,  the  code 
calculates  coefficients  for  analytical  functions  to  fit 
the  Monte  Carlo  data.  The  purpose  of  the  functional  repre- 
sentation of  the  data  is  twofold:  first,  to  provide  useful 

expressions  to  describe  the  profiles,  and  second,  to  smooth 
out  statistical  fluctuations  in  the  Monte  Carlo  data.  The 
computed  charge  and  energy  current  profiles  are  fitted  to 
functions  of  the  form 

y(x)  = A exp(Bx  + Cx2  + Dx3)  (8) 

The  expressions  for  the  charge  and  energy  deposition  pro- 
files are  determined  through  analytical  differentiation  of 
the  functions. 
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2.3  SPECIAL  PURPOSE  CODES 

In  addition  to  the  emission  code  and  interface  code 
described  above,  several  other  special  purpose  POEM  codes 
have  been  developed.  These  codes  are  basically  the  same 
as  the  previous  codes,  generally  differing  only  in  the 
random  walk  and  result  subroutines.  The  applications  of 
these  codes  are  described  briefly  in  the  following  paragraphs. 
The  user  instructions  for  these  codes  are  not  included  in 
this  report,  nor  are  the  appropriate  subroutines  included 
in  the  released  POEM  library.  The  instructions  and  rou- 
tines, however,  may  be  obtained  from  the  author. 

2.3.1  Multislab  Transition  Zone  Dose  and  Current:  (MS) 

This  code  calculates  the  photo-Compton  current  and 
dose  profiles  through  a stack  of  thin  foils,  each  foil 
thinner  than  the  maximum  electron  range.  The  problem  geo- 
metry is  defined  in  Figure  4.  A stack  of  up  to  20  planar 
slabs  of  up  to  four  different  materials  is  exposed  to  plane 
wave  photon  radiation  incident  at  an  angle  6q  with  respect 
to  normal  to  the  slabs.  The  configuration  is  zoned  into 
up  to  100  deposition  zones.  The  photo-Compton  charge  and 
energy  current  and  deposition  profiles  are  calculated. 

Sample  results  are  included  in  Reference  14. 

2.3.2  Scanning  Electron  Microscope  Dose  and  Current:  (SEM) 

This  code  calculates  the  two-dimensional  current  and 
dose  profiles  in  a one  to  three  layer  structure  irradiated 
by  the  electron  beam  of  a scanning  electron  microscope  (SEM) . 
The  problem  geometry  is  defined  in  Figure  5.  A structure 
composed  of  up  to  three  thin,  planar  slabs  of  material  is 
normally  irradiated  by  a narrow  electron  beam.  The  code 
uses  a cylindrical  coordinate  system  with  the  z-axis  along 
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Figure  5 


PROBLEM  GEOMETRY:  SEM  CURRENT  AND  DOSE 


the  electron  beam.  The  code  calculates  the  two-dimensional 
current  and  dose  distributions,  j(r,z)  and  D(r,z),  the  back- 
scattered  charge  and  energy,  and  the  radial  distribution  of 
the  backscattered  electrons.  A discussion  of  the  code  and 
sample  results  for  MOS  devices  are  included  in  Reference  15. 

2.3.3  Transition  Zone  Dose  and  Current  in  Dielectrics: 

(DIED 

This  code  calculaces  the  transition  zone  dose  and 
current  in  an  x-  or  yirradiated  dielectric  near  a planar 
interface.  The  electron  transport  calculation  includes  the 
electron  acceleration  in  the  electric  field  of  the  trapped 
charge  in  the  diexectric.  The  time-dependent  field  calcu- 
lation includes  a time-dependent , spatially-dependent 
treatment  of  the  transient  conduction  currents  in  the 
dielectric.  The  problem  geometry  is  shown  in  Figure  6. 

A planar  interface  between  a conductor  and  dielectric  is 
exposed  to  plane  wave  x or  y radiation.  The  xy-plane  lie^ 
along  the  interface;  the  z-axis  is  directed  into  the  dielec- 
tric. The  code  calculates  the  photo-Compton  current  j(z,t), 
the  true  current  J(z,t), 

J(z,t)  = j(z,t)  + o ( z , t ) E ( z , t ) (9) 

the  charge  density 

t 


o 

the  radiation-induced  conductivity  o(z,t),  and  the  electric 
field  E(z,t) . 
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PROBLEM  GEOMETRY:  TRANSITION  ZONE  DOSE 

AND  CURRENT  IN  DIELECTRICS 


A discussion  of  the  dielectric  code  and  sample  re- 
sults are  found  in  Reference  16.  The  code  is  presently 
being  modified  to  treat  a planar  conductor/dielectric/con- 
ductor sandwich  configuration.  It  is  planned  to  build  a 
three-dimensional  version  of  this  code  within  the  next 
year . 

2.3.4  I.EMP  in  a Plated  Wire  Memory:  (PWM) 

A three-dimensional  version  of  POEM  was  developed 
for  the  calculation  of  the  induced  bit  wire  and  word  strap 
currents  in  a plated  wire  memory  exposed  to  x or  y radiation. 
The  electron  transport  is  calculated  for  the  problem  geometry 
defined  in  Figure  7.  The  charge  displacement  in  the  mem- 
ory structure  is  translated  into  induced  currents  on  the 
bit  wires  and  word  straps  using  a Green's  function  compu- 
tational technique.  The  coupling  is  represented  by  equiva- 
lent current  sources  driving  a lumped  element  circuit  repre- 
sentation of  the  memory  electrical  configuration. 


2.3.5  Compton  Currents  in  a Right  Circular  Cylinder 


A two-dimensional  version  of  POEM  was  developed  for 
the  calculation  of  the  Compton  currents  and  gas  ionization 
in  a gas-filled  right  circular  cylinder  irradiated  end-on 
by  photons  (Figure  8)  . The  two-dimensional  distributions 
of  the  axial  and  radial  space  current  components  and  the 
energy  deposition  in  the  gas  are  calculated  for  an  arbitrary 
incident  photon  energy-angular  distribution.  The  electron 
sources  include  the  front,  back  and  side  walls  of  the  cylin- 
der, the  gas,  and  an  optional,  coaxial  inner  cylinder. 
Space-charge  limitation  is  negligible  for  the  application 
for  which  the  code  was  designed,  so  self-consistent  effects 
are  ignored  in  the  transport  calculation.  These  transport 


I 

I 


Figure  7 

PROBLEM  GEOMETRY  - PLATED  WIRE  MEMORY 
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calculations  are  used  as  inputs  to  a source  region  EMP  code 
for  the  calculation  of  the  fields  within  the  cavity  and  the 
surface  current  induced  on  the  inner  cylinder.  The  cods  and 
sample  calculations  are  discussed  in  Reference  29. 
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The  computation  of  the  x-ray  photoemission  or  the 
transition  zone  dose  and  current  involves  three  main  parts: 

(1)  calculation  of  the  photon-generated  electron  source 
distribution,  (2)  calculation  of  the  electron  transport, 
and  (3)  processing  of  the  transport  data  to  obtain  the 
results.  In  the  following  paragraphs  the  computational 
procedures  used  in  each  of  these  calculations  are  discussed. 
In  order  to  efficiently  provide  an  overview  of  these  pro- 
cedures, detailed  discussion  of  the  cross  sections  employed 
in  the  calculations  is  postponed  until  Section  4. 

3.1  ELECTRON  SOURCE  CALCULATION 

The  calculation  of  the  photon-generated  electron 
source  distribution  is  performed  using  an  analytical  model. 

To  do  this,  it  is  assumed  that  (1)  the  photon  flux  (intensity 
and  spectrum)  is  defined  at  the  emission  surface  (or  material 
interface) , and  (2)  the  photon  flux  is  constant  over  the 
region  of  the  material  within  an  electron  range  from  the 
surface.  This  means  that  it  is  assumed  that  the  maximum 
electron  range  is  small  compared  to  the  photon  mean-free 
path, 


r 
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1 

u 


a good  assumption  for  the  range  of  applicability  of  the  code. 
In  this  approximation  the  electron  source  density  is  constant 


over  the  region  contributing  to  the  electron  emission.  We 
call  this  region  the  source  region. 

Let  the  thickness  of  the  electron  source  region  be 
hz  (g/cm2 ) (Az<re) , and  the  planar  fluence  of  the  photon 
group  of  energy  hv  by  y(hv) (photons/cm2).  For  a Compton 
cross  section  o^nc(cm2/g),  the  density  Sc (electron/cm2 ) of 
Compton  electrons  produced  in  the  source  region  by  photons 
of  energy  hv  is 

SC(hv)  = Y°inc  Az  (11) 

For  a photoelectric  cross  section  a , the  density 
of  pnotoelectrons  is 

(hv)  = yotAz  (12) 

For  each  photon  energy  group,  the  source  calculation  includes 
two  photoelectron  energy  groups  as  defined  in  Table  3 where 
Ek  is  the  K photoelectric  absorbtion  edge  and  E^-,  E^,  and 
E—  are  the  mean  L,  M,  and  N edges. 

The  mean  edges  and  the  relative  probabilities  for 
that  shell  in  which  the  ionization  occurs  have  been  calcu- 
lated using  the  ratios  of  the  photoelectric  cross  sections 
just  below  and  just  above  the  absorbtion  edge.  The  photo- 
electric cross  section  undergoes  a jump  at  each  atomic 
shell  binding  energy  E^  (edge) . The  ratio  of  the  cross 
section  just  below  the  edge  E^  to  the  cross  section  just 
above  the  edge  is  the  probability  that  the  ionization  occurs 

in  one  of  the  shells  with  binding  energy  less  than  E,  . We 
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make  the  common  assumption  that  this  probability  is  con- 
stant for  all  photon  energies  between  E^  and  the  next  larger 
shell  binding  energy.  For  example,  for  hv>E  the  probability 

— Jt\ 
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that  the  ionization  occurs  in  the  L shell  or  lower  binding 
energy  shell  is  rK  where  rK  is  the  ratio  of  the  cross  sec- 
tions just  above  and  just  below  the  K edge.  The  probability 
that  the  ionization  occurs  in  the  K shell  is  1-r^.  We  use 
only  two  photoelectron  groups;  we  lump  the  L and  lower 
binding  energy  shell  ionizations  together  by  defining  the 
relative  K and  L shell  ionization  probabilities  to  be 

p = 1 — r (13) 

PK  K 


P 


L 


(14) 


for  hv>EK. 


Similarly  for  E— <hv<E^  we  define 


PM  = 1 rLirL2rL3 


PIu  1 rLi 


PL2  rLi (1  rL2} 


L2 


(1-r, 
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where  rT  ■ is  the  cross  section  ratio  at  the  L,  edge.  We 
LI 

lump  the  three  L subshell  ionizations  together  by  defining 
a relative  L shell  ionization  probability, 


= PL i +PL  2+PL  3 


rLirL2rL, 


(19) 
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and  a mean  L edge, 


3 


i=l 


Li  Li 


(20) 


The  mean  I,  edge  is  used  in  calculation  of  the  L photoelectron 
kinetic  energy  as  shown  in  Table  3.  Similarly,  for 
E— <hv<E—  we  define  relative  ionization  probabilities  pM  and 
PN  and  a mean  M edge  Ej^: 

5 

Pm  PMi  (21> 

i=l 
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i=l 


(22) 


where 
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Mi-1  Mi 


(23) 


The  photoelectron  source  group  definition  is  summar- 
ized in  Table  3.  The  mean  N edge  is  approximated  by  zero. 
The  code  is  invalid  for  hv<E-  (Ej^<5  keV)  . The  relative 
ionization  probabilities  defined  in  Table  3 and  the  mean 
edges  defined  in  Equations  (20)  and  (22)  have  been  precalcu- 
lated and  stored  on  the  photon  interaction  data  tape. 


Following  the  photoelectric  ionization  the  atom  will 
relax  emitting  either  a fluorescent  photon  or  an  Auger 
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electron.  For  the  electron  source  calculation,  one  Auger 
electron  energy  group  is  included  for  each  photon  energy 
group  as  defined  in  Table  4.  The  K shell  fluorescence 

_ 19 

yield  <*,.  have  been  obtained  from  the  Storm  and  Israel  and 
L 25 

Bambynek,  et  al.  compilations,  with  preference  given  to 
the  Bambynek  compilation,  and  stored  on  the  photon  inter- 
action data  tape.  The  mean  M fluorescence  yield  is 
approximated  by  zero. 

The  total  electron  source  density  S'(hv)  produced 
by  the  photon  group  of  energy  hv  is  the  sum  of  the  Compton, 
photoelectric,  and  Auger  electrons, 

S' (hv)  = Sc(hv)  + ST(hv)  + SA(hv)  (24) 


The  total  source  density  S is  obtained  by  summing  S' (hv) 
over  the  photon  spectrum: 


120 

S = ^ S'  (hvjL) 
i=l 


(25) 


The  photon  spectrum  is  defined  in  the  code  as  (1)  a mono- 
chromatic spectrum,  (2)  a black  body  spectrum  internally 
calculated  over  120  photon  energy  groups,  or  (3)  an  arbi- 
trary spectrum  defined  over  up  to  120  photon  energy  groups 
The  spectrum  is  internally  normalized  to  a planar  fluence 
of  1 photon/cm2. 

The  energy-angular  distributions  of  the  Compton, 
photo-,  and  Auger  electrons  are  calculated  using  the  cross 
sections  discussed  in  Section  4.  These  distributions  are 
summed  over  the  photon  spectrum  to  obtain  a joint  energy- 
angular  distribution  function  for  the  electron  source. 
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This  is  a discrete  distribution  function  defined  over  a 

two-dimensional  array  of  20  energies  by  20  angles.  The 

function  is  stored  as  an  array  H(E.,6.),  i=l,20,  j=l,20 

3 

where 


H(E.,0.)  = p (E i ) + p(E,)  + ...  + P(Ei_1) 


+ p(Ei) [p(0, |E±)  + p(02 |E.) 


+ ...  + p (0 j | E^)  ] 


(26) 
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The  angular  distribution  is  calculated  with  respect  to  the 
photon  direction  of  incidence  (Figure  9) . The  polar  an- 
gular (0)  distribution  is  defined  ov°r  0 to  u.  The  azimu- 
thal angular  (4>)  distribution  about  the  direction  of 
photon  incidence  is  assumed  to  be  uniform  over  0 to  2tt . 

In  addition  to  the  joint,  cumulative  distribution  function, 
the  marginal  (one-dimensional)  energy  and  angular  distri- 
bution functions  p(E^)  and  ptO^),  i-i,20,  are  calculated 
for  the  output. 
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Figure  9 


COORDINATE  SYSTEM  x'y'z'  FOR  ELECTRON  SOURCE 
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The  electron  source  density  S and  the  energy-angular 
distribution  H(E^,0j)  define  the  electron  source  for  input 
to  the  electron  transport  calculation. 

3.2  ELECTRON  TRANSPORT 

POEM  calculates  the  electron  transport  using  the 
Monte  Carlo  method.  The  transport  is  simulated  by  tracing 
the  trajectories  of  many  particles  through  a random  walk 
calculation.  The  initial  conditions  for  each  particle  - 
energy,  direction,  and  location  - are  random  selected  from 
the  electron  source  probability  distributions.  The  weight 
of  each  particle  is  determined  by  the  electron  source 
density.  The  trajectory  of  the  particle  is  traced  using 
precalculated  slowing  down  and  scattering  distributions. 
Relevant  quantities  - energy  and  charge  currents  at  zone 
boundaries  and  energy  and  charge  depositions  within  zones  - 
are  scored  for  each  trajectory.  Estimates  of  the  desired 
output  quantities  and  estimates  of  the  statistical  errors 
are  made  by  averaging  the  results  for  many  particle 
histories . 

3.2.1  The  Monte  Carlo  Method 

As  an  example  of  the  Monte  Carlo  method  used  in 
POEM,  we  discuss  the  calculation  of  the  electron  emission 
from  a planar  surface.  Extension  of  the  method  to  the 
transition  zone  dose  and  current  calculation  is  straight- 
forward. We  express  the  yield  in  electrons/photon-Sr-keV 

as  a function  y 1 (E . , 0 . , <J>,  ) over  a three-dimensional  array 

3 ^ 

of  n.  energy  bins,  n.  6-angle  bins,  and  n,  <{>-angle  bins. 

1 3 k 

Each  bin  is  characterized  by  its  midpoint,  for  example, 
the  energy  bin  E^  contains  the  energy  interval  (E^  - AE/2, 

E^  + AE/2).  Let  Ne  (electrons/cm2)  be  the  density  of 
electrons  produced  in  the  source  region.  Let  F(E,0,(|>,z) 
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be  the  cumulative  probability  distribution  function  for 
the  initial  energy  E,  direction  (9  ,<{>),  and  location  z of 
the  electrons.  The  yield  is  then  given  by 


y'  (Ei'°j'V  = NeyP(Ei,9j,<j>k|E,0><J,,z) 


dF  (E  , 0 , <(» , z) 


(27) 


where  p (E^  , 9 ^ J E , 9 , <J> , z)  is  the  probability  that  an  elect- 
ron with  initial  energy  E,  direction  (9,<J>),  and  location 
z is  emitted  into  energy  bin  E^,  9-angle  bin  9^,  and  <j>- 
angle  bin  <})^.  Let  the  vector  X = (E,0,4>,z)  define  the 
initial  electron  conditions.  We  denote  y„k  = y 1 (E^ , 9 ^ , 4>^) 
and  p.  (X)  = p (E  . , 9 . , <f>,  I E,  6 , <f> , z)  . T’he  equation  can  be 

1]K  ~ 1 J K 

rewritten  as 


*ijk  * Ke/ Pijk(X)dF(X) 


(28) 


POEM  uses  the  Monte  Carlo  method  to  evaluate  the 
integral  in  Equation  (28)  . Particles  are  started  with  ini- 
tial conditions  selected  from  the  source  distribution  F(X), 
and  the  particle  trajectories  are  traced  until  the  particle 
is  emitted  from  the  surface,  or  until  the  residual  range 
of  the  particle  is  insufficient  for  the  particle  to  reach 
the  surface.  For  particle  number  n we  define  the  random 
variable  Cijk(n)  by 


1,  if  particle  n is  emitted 

with  conditions  E.,9.,4>, 

JL3^ 

0,  otherwise 


(29) 
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We  define  the  random  variable  Y. (N)  as 

13k 


vi]k(N)  ■ 


N 

= irX)  ci:*<n 


) . 


(30) 


n=  _ 


We  calculate  Y...  (N)  by  tracing  the  trajectories  of  N 

1 J K 

particles.  Each  time  a particle  is  emitted  with  conditions 
(E  . , b . , <J>.  ) , we  increment  Y...  (N)  by  N /N.  The  expected 

-i-jK  1 ] K G 

value  of  Y.  (N)  is  y . . , : 

13k  -*13^ 


I I Ne 

E Y.  ..  (N)  = — 

I 13kv  I N 


N 


E 

n=l 


E ^ijk<n> 


e 

N 


'f  pii' 


N / d.  (X)  dF  (X) 


= y 


ijk 


(31) 


The  calculated  quantity  Y^^N)  is  called  an  unbiased  esti- 
mator of  the  yield  The  law  of  large  numbers  tells  us 
that  for  any  e>o, 


Lirr^  prob||Yijk(N)  - yijkl<e|  = 1 


(32) 


That  is,  Yijk(N)  tends  to  y^jk  in  probability.  Thus,  if 
we  calculate  a sufficient  number  of  particle  histories, 
we  can  estimate  the  yield  y^k  by  the  estimator  Y^k(N)  to 
to  within  arbitrary  precision. 

The  variance  of  the  estimator  Y. (N)  is 

13k 


VAr|y.  (N)  ! = E J(Y.  (N)  - y.  ) 2<  = J 
I 13k  ) r 13k'  j 13k7  n p 


(33) 
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where 


3p ‘/tpijk(s)  ‘ yijk'2dF(s>  ,34) 

is  the  variance  of  the  emission  probability  distribution 
p...  (X).  The  square  root  of  the  variance  of  Y.  ,,(N)  is  the 

ljK  — 13K 

standard  error 


Oy(N) 


(35) 


Equation  (35)  shows  that  the  statistical  error  a (N)  is  in- 
versely proportional  to  the  square  root  of  the  number  of 
particle  histories.  Thus,  to  reduce  the  error  by  a factor 
k requires  increasing  the  number  of  particle  histories  by 
a factor  of  k2. 


The  standard  error  oy  is  not  a computed  quantity  in 
the  Monte  Carlo  calculation.  However,  just  as  we  obtain  an 
estimator  for  the  yield  y. ..  we  obtain  an  estimator  for  the 

13K 


standard  error, 
s4  is  given  by 


The  unbiased  estimator  for  is  s where 


s2  = 


1 

N-l 


IN 

n=l 


<T(X  > - Y?.,  (N) 


— n 


i jk 


(36) 


3.2.2  Variance  Reduction 

Using  the  Monte  Carlo  method  we  can  obtain  predic- 
tions for  the  yield  y 1 (E . , 0 ■ , ) and  the  statistical  error 

X J K 

of  the  estimated  yield.  The  problem  is  that  emission  of 

a particle  into  any  one  of  the  bins  (E.,9.,<J>,)  of  the  three- 

1 3 x 

dimensional  emission  distribution  is  an  infrequent  event. 

For  example,  for  back  emission  from  aluminum  irradiated  by 
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10  keV  x-rays,  the  mean  probability  of  emission  of  a source 
region  electron  is  0.058.  The  mean  probability  for  the 
emission  of  the  electron  into  one  of  20  bins  of  a one- 
dimensional angular  distribution  is  0.0029.  The  mean  prob- 
ability for  emission  into  one  of  the  160  bins  of  a two- 
dimensional  angular  distribution  over  20  0-bins  by  8 <j>-bins 
is  0 00036.  This  means  that  in  a straightforward  Monte 
Carlo  calculation  one  would  expect  on  the  average  four 
particles  out  of  10,000  histories  to  be  emitted  into  any 
one  bin  of  the  two-dimensional  angular  distribution.  To 
obtain  a statistically  meaningful  distribution  would  re- 
quire many  hundred- thousand  histories  at  considerable  com- 
putation expense. 

Since  to  halve  the  standard  error  one  must  quadruple 

the  number  of  particle  histories,  error  reduction  through 

increased  numbers  of  particle  histories  is  inefficient.  A 

more  efficient  method  of  error  reduction  is  to  reduce  the 

standard  deviation  a of  the  parent  distribution  rather  than 

P 

increase  the  number  of  histories  N [Equation  (35)].  This 
can  be  done  through  variance  reduction  techniques. 

Variance  reduction  in  Monte  Carlo  can  be  viewed  as 
the  introduction  of  known,  usually  qualitative,  information 
about  the  process  being  simulated  into  the  Monte  Carlo  simu- 
lation in  an  explicit  and  quantitative  manner.  Through  the 
introduction  of  such  information,  the  amount  of  information 
which  must  be  obtained  through  the  simulation  is  reduced. 
Variance  reduction  techniques  for  radiation  transport  cal- 
culations are  highly  developed  - through  judicious  applica- 
tion of  these  techniques  code  efficiencies  can  in  some 
cases  be  increased  by  factors  of  up  to  several  thousand. 

A good  mathematical  treatment  of  variance  reduction  is 
found  in  Reference  17;  a convenient,  practical  discussion 
is  found  in  Reference  18. 
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The  variance  reduction  technique  which  has  been  used 
to  most  advantage  in  POEM  is  importance  sampling.  The  elect- 
rons which  contribute  most  significantly  to  the  emission  are 
sampled  preferentially  through  biasing  the  electron  source 
distributions.  The  weights  of  the  transport  particles  are 
adjusted  accordingly  so  as  not  to  bias  the  estimator.  Im- 
portance sampling  increases  the  number  of  particle  emission.*? 
and  decreases  the  computation  time  wasted  on  those  particles 
which  do  not  contribute  to  the  emission  distribution. 

Through  importance  sampling  the  POEM  code  efficiency  has 
been  increased  by  about  a factor  of  20;  i.e.,  to  obtain  the 
same  standard  error  the  code  with  importance  sampling  re- 
quires 1/20  the  computation  time  required  for  straightforward 
Monte  Carlo. 

To  illustrate  the  importance  sampling  concept,  con- 
sider a one-dimensional  simplification  of  the  integral  in 
Equation  128) : 


y = n 


p (x)  dF (x) . 


Where  y is  yield,  p(x)  is  the  probability  of  emission  given 
initial  condition  x and  F(x)  is  the  cumulative  probability 
distribution  for  x.  Rather  than  compute  this  integral,  we 
compute  the  integral 


y = N 


■<*> 


where 


dF'  (x)  = g(x)dF(x)  = g(x)f(x)dx 
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I 
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Ve  track  N particles  selecting  their  initial  condition  x 
from  the  cumulative  probability  distribution  function  P' 
rather  than  F and  weight  the  particles  ry  l/g(xn).  Just 
as  before  we  define  an  estimator  Y'(N): 

N 

N m 1 

Y 1 (N)  = / W £,  (4 

N / j n ^n 

n=l 


where 


n 


1,  if  particle  n is  emitted 
0,  otherwise 


(41) 


and 


f. 

I 

t 


t 


t 

( 


W 


n 


g(xn) 


The  expected  value  of  the  estimator  is 


EY'(N)  = 


N 

e 

N 


N 


E 

n=l 


Wn1 


N 

e 

N 


N/ 


P (x) 
g(x) 


dF' (x) 


= y 


(42) 


Thus,  Y ' (N)  is  an  unbiased  estimator  of  y.  The  standard 
error  Oy , is 


a 


Y' 


1 

M °p/g 


(43) 


61 
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where  °p/g  is  the  standard  deviation  of  distribution  p(x)/g(x). 


If  g(x)  is  a reasonable  approximation  of  p(x)  over  the  range 

of  x,  i.e.,  g(x)p(x),  x . -'x-x  , then  the  standard  devi- 

3 r nun-  --  max 

will  be  small  compared  with  the  standard  devia- 


ation 


p/g 


tion  of  the  distribution  p(x).  Therefore,  comparing 


Equations  (35)  and  (43),  the  standard  error  for  the 

calculation  using  importance  sampling  will  be  small  compared 


with  the  standard  error  oy  for  the  calculation  using  the 


straightforward  Monte  Carlo. 

To  incorporate  importance  sampling  into  the  emission 


calculation  we  bias  the  sampling  distribution  f (x)  for  the 
electron  initial  conditions  with  a functior  g(x)  which  is 
the  best  approximation  obtainable  o:  the  probability  of 
emission  p(x).  In  POEM  importance  sampling  is  used  in  the 
selection  of  the  particle  initial  location  and  initial  energy. 

For  electrons  of  initial  energy  E and  range  R(E), 
emission  is  possible  only  for  those  electrons  with  initial 
location  at  distance  from  the  surface  2-"R(E).  Thus,  in 
POEM,  initial  location  is  sampled  only  over  the  range 
0^z<R(E).  The  particles  are  weighted  by  W = R(E)Az  where 
Lz  is  the  thickness  of  the  source  region.  Because  of 
scattering,  electrons  with  a given  range  R(E)  which  start 
nearer  the  surface  have  higher  probability  of  emission.  In 
POEM,  the  sampling  distribution  for  initial  electron  distance 
z from  the  surface  is  biased  by  a linear  function  g(z)  which 
is  zero  at  z = R(E)  and  maximum  at  z = 0.  The  particle 
weight  is  multiplied  by  a factor  l/g(z). 

Within  the  source  region,  electrons  with  higher 
initial  energy  have  higher  probability  of  emission.  In 
POEM,  the  sampling  distribution  for  initial  electron  energy 
is  biased  towards  high  energy.  The  biasing  and  weighting 
are  performed  in  the  random  selection  procedure  as  explained 
in  Section  3.2.4. 
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Importance  sampling  is  not  used  in  POEM  for  the 
selection  of  initial  electron  direction.  Clearly  impor- 
tance sampling  could  be  advantageous  here  - those  electrons 
starting  in  a direction  towards  the  surface  have  higher 
probability  of  emission  than  those  starting  in  a direction 
away  from  the  surface.  Appropriate  biasing  of  the  sampling 
distribution  for  initial  direction  is  currently  being  in- 
vestigated . 

Another  form  of  importance  sampling  used  in  POEM 
is  particle  killing.  A particle  whose  residual  range  is 
less  than  the  minimum  distance  to  an  emission  surface  (or 
deposition  zone  boundary)  is  killed  since  it  cannot  contri- 
bute to  the  distribution.  This  is  an  obvious  variance  re- 
duction technique  which  is  also  available  in  such  codes  as 
ETRAN  and  SANDYL. 

Other  variance  reduction  techniques  are  being  tried 
in  experimental  versions  of  POEM.  These  include  track 
splitting,  control  variates,  and  orthonormal  polynomial 
expansions  of  output  distributions.  These  techniques,  if 
proven  effective,  will  be  incorporated  into  POEM  and  reported 
in  subsequent  publications. 

3.2.3  Transport  Grid 

POEM  is  a multiple-scatter  electron  transport  code. 
Rather  than  directly  simulate  each  collision,  as  is  commonly 
done  in  photon  transport  codes,  each  step  of  the  random  walk 
calculation  simulates  the  effects  of  many  collisions.  The 
step  sizes,  the  energy  loss  increments,  and  the  multiple 
scattering  probability  distributions  corresponding  to  each 
step  of  the  random  walk  are  predetermined  using  the  code 
DATAGEN . 
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An  energy  grid  is  set  up  using  logarithmic  energy 
steps  from  the  maximum  electron  energy  down  to  1 keV. 
Typically,  the  energy  step  sizes  are  set  such  that  an  elect- 
ron loses  one-half  of  its  energy  in  eight  steps.  The  number 
of  steps  to  halve  tbe  energy  is  NCYC.  The  energy  of  each 
successive  grid  point  is  reduced  by  a factor  EFAC  where 
EFAC**NCYC  = 1/2: 

E(N)  = E(N-l)  * EFAC.  (44) 

At  low  energies,  E-_  10  keV,  it  is  found  that  the  energy  steps 
are  needlessly  small  for  NCYC  = 8,  particularly  for  high-Z 
materials.  The  corresponding  spatial  step  sizes  are  so 
small  that  there  is  a sizable  probability  of  zero  deflection 
along  the  step.  Consequently,  the  logarithmic  spacing  of 
the  energy  grid  is  modified  below  about  10  keV.  Typically 
we  set  NCYC  = 8 at  high  energies,  change  to  4 at  about 
10  keV,  and  to  2 at  about  5 keV. 

To  insure  sufficient  detail  in  the  random  walk  cal- 
culation, the  major  energy  grid  is  subdivided  into  NSUB 
substeps : 


AE  (N) 


E (N)  - E (N+l ) 
NSUB 


(45) 


The  appropriate  value  for  NSUB  depends  on  the  transport 
medium  atomic  number,  and  also  to  some  extent  on  the  user 
application.  Higher  atomic  number  materials  require  more 
substeps.  We  have  obtained  good  results  for  NSUB  values 
typified  by  Table  5.  (These  values  are  somewhat  smaller 
than  those  given  by  the  empirical  formula  recommended  for 
SANDYL.2) 
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Table  5 

TYPICAL  NSUB  VALUES 


MATERIAL 

Z 

NSUB 

C 

6 

2 

Al 

13 

4 

Sr. 

50 

6 

Au 

79 

8 

POEM  uses  the  continuous  slowing  down  (csda)  approxi- 
mation: the  spatial  step  size  A s is  uniquely  determined  by 

the  energy  step  size  AE.  (POEM  neglects  energy  loss  strag- 
gling. This  does  not  appear,  however,  to  be  an  important 
effect  in  the  x-ray  photoemission  or  transition  zone  dose 
applications  for  which  POEM  was  designed  - we  obtain  good 
agreement  with  experiment.  Thus,  the  inclusion  of  straggl- 
ing does  not  at  this  time  appear  to  warrant  the  cost,  i.e., 
the  attendant  degradation  of  code  efficiency.)  The  spatial 
step  size  As  is  calculated  from  the  slowing  down  distribution 
dE/ds (E) : 


E-AE 


dE 

-dE/ds 


The  stopping  power  -dE/ds  includes  both  collision  and  radi- 
ation terms 


dE  _ /dE  \ / dE  \ 

ds  ' OOcol  rOrad' 


I 


In  traversing  the  pathlength  increment  As,  the  elec- 
tron undergoes  a number  of  single  scatterings.  The  proba- 
bility distribution  for  the  total  deflection  of  the  electron 
is  calculated  using  the  Goudsmi tt-Saunderson  multiple 
scattering  theory.  One  distribution  is  calculated  for  each 
point  N of  the  major  energy  grid.  A cumulative  probability 
distribution  G^(b)  is  calculated  for  an  electron  of  energy 
E^  for  energy  loss  AE  and  pathlength  increment  As.  The 
distribution  is  calculated  over  an  angular  grid  defined  by 
the  endpoints  6 = 0.5,  1.0,  1.5,  2,  3,  4,  6,  8,  10,  15,  20, 
30,  40,  60,  80,  100,  120,  140,  160,  and  180  degrees. 

The  transport  grid  - the  major  energy  grid  endpoints 
E^,  the  energy  substeps  AE^,  the  pathlength  increments  As^, 
and  the  multiple  scattering  distributions  GN(0)  are  calcu- 
lated by  the  code  DATAGEN  and  are  stored  for  use  in  the 
random  walk  calculation. 


3.2.4 


Random  Walk 


The  random  walk  calculation  for  each  particle  starts 
with  random  selection  of  the  initial  energy,  direction,  and 
location  and  proceeds  through  spatial  steps  and  angular 
deflections  until  the  particle  escapes  the  transport  region 
or  is  absorbed.  Along  the  particle  trajectory,  zone  and 
material  boundary  crossings  are  scored.  In  the  following 
paragraphs,  we  discuss  the  procedures  for  trajectory  ori- 
gination, trajectory  construction,  trajectory  termination, 
boundary  crossings,  and  scoring. 

Trajectory  Origination  - The  electron  history  starts 
with  random  selection  of  the  initial  energy  EQ.  The  cumu- 
lative distribution  function  F„(E)  is  determined  from  the 
cumulative,  two-dimensional  source  distribution  H(E,6) 
(Equation  (26) ] : 
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Fe(E)  = H (E , Tt ) . 


(48) 


The  energy  Eq  is  selected  by  generating  a random  number  rQ 
between  0 and  1 and  setting  E = F;1 (r  ) where  F ~ • is  the 
inverse  of  the  distribution  function  F . The  procedure  is 

Ij 

illustrated  in  Figure  10.  If  rQ  were  selected  from  the 
uniform  distribution  over  the  interval  0 to  1,  U(0,1),  then 
the  selected  distribution  of  initial  energies  would  simu- 
late the  source  distribution;  i.e.,  after  many  histories, 
the  selected  initial  energy  cumulative  distribution  would 
approach  F^  (E)  . We,  however,  want  to  bias  the  sampling 
distribution  towards  high  energy.  We  select  a random  num- 
ber u from  U(0,1)  and  then  set  r = 1 - u2 . This  bias 
o o o 

was  found  by  experience  to  maximize  one  code  efficiency 
over  the  range  of  application.  The  appropriate  particle 
weight  factor  can  be  shown  to  be 

WE  = 2uq . (49) 

Note  that  the  initial  energy  is  not  selected  from  a discrete 
distribution  but  from  a continuous  distribution.  This  is 
done  by  linear  interpolation  in  the  stored,  discrete  source 
distribution,  as  is  indicated  graphically  in  Figure  10. 

Next,  the  initial  distance  zQ  from  the  surface  is 
selected.  The  source  density  is  assumed  uniform  over  the 
source  region,  0<_zq<Az.  The  electror  of  initial  energy  EQ 
can  only  be  emitted,  however,  from  the  region  0<z  <_R(EQ). 
Thus,  the  sampling  is  restricted  to  the  region  0£zq_<R(Eo)  . 
Accordingly,  the  particle  weight  is  modified  by  the  factor 
R(Eq)/Az.  We  want  to  bias  the  sampling  distribution  towards 
zQ  - 0 . A linear  sampling  distribution  is  used  which  is 
maximum  at  z = 0 and  zero  at  zQ  = R(F.).  This  is  effected 
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RANDOM  SELECTION  OF  INITIAL  ENERGY 
FROM  CUMULATIVE  PROBABILITY  DISTRIBUTION 


by  selecting  a random  number  uq  from  U(0,1),  setting  rQ  = 
v'u^,  and  setting 

zQ  = (1  - rQ)R(Eo) . (49) 


The  appropriate  weighting  factor  is  l/2rQ.  The  total 
weighting  factor  Wz  for  the  z selection  is  then 


W 


z 


RUy 

Lz 


(59) 


Selection  of  the  initial  energy  Eq  determines  the 
distribution  function  for  the  initial  6-angle.  Let  the 
selected  initial  energy  Eq  be  in  the  energy  group  EN  of  the 
cumulative  source  distribution  function  H(E,6).  The  cumu- 
lative distribution  function  for  the  initial  6-angle 
(with  respect  to  the  photon  direction  of  incidence)  is 


H(E  ,6)  - H(E  ,0) 
Fe(0)  = II(En,tt)  - H(EN,0)  ' 


(51) 


No  biasing  is  used  in  the  selection  of  the  angle  0Q.  A ran- 
dom number  uq  is  selected  from  U(0,1).  The  initial  angle 
, with  respect  to  the  photon  direction  of  incidence,  is 
selected  using  the  same  procedure  as  for  the  selection  of 
the  initial  energy  EQ: 


O'  just  as  E , is  selected  from  a continuous  distribution 
o J c 

through  linear  interpolation  in  the  computed,  discrete 
source  distribution. 
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An  initial  azimuthal  angle  (j>^  is  selected  from  a 

uniform  distribution  about  the  photon  direction  of  incidence. 

A random  number  u is  selected  from  U(0,1),  and  <(>'  is  set 

o 0 

<\ ' = u • 2ir  . (53) 

o o 

The  initial  direction  ((V,  <J>^)  is  selected  from  angu- 
lar distributions  defined  with  respect  to  the  photon  direc- 
tion of  incidence.  Transformation  to  the  direction  (0O/4>O) 
in  the  problem  geometry  coordinate  system  is  performed  using 
a straightforward  orthogonal  transformation. 

The  total  weight  Wq  of  the  particle  is  the  product 
of  the  weight  factors  W£  and  W^, 

w = w„w  . (54) 

O E Z 


Trajectory  Construction  - A particle  of  weight  WQ  and 
initial  energy  Eq  starts  out  at  location  zQ  and  direction 
(0  , 0 ).  The  particle  takes  a partial  step  to  reduce  its 
energy  to  the  nearest  lower  point  of  the  transport  energy 
grid.  Let  the  initial  energy  Eq  be  bound  by  the  major 

energy  grid  points  E^  and  EN+^,  ^N— Eo>^N+l'  an<^ 
the  energy  loss  increment  required  to  reach  the  next  sub- 
division point  of  the  energy  grid.  Then  the  partial  step 
size  6s  is  calculated  by  linear  interpolation: 


6s 


AsN  ’ 


(55) 


Thereafter,  the  calculation  steps  through  the  trans- 
port grid,  incrementing  the  pathlength  and  decreasing  the 
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| particle  energy.  At  the  end  of  each  pathlength  increment 

l As  , t^ie  c°de  samples  the  angular  deflection  AO  , with  re- 

^ i*  IN 

l spect  to  the  present  direction  of  travel,  from  the  multiple 

I - scattering  distribution  G^(0).  The  angle  is  calculated  by 

selecting  a random  number  uq  from  U(0,1)  and  setting 


AO 


N 


(V 


(56) 


The  angular  deflection  is  sampled  from  a continuous  distri- 
bution through  linear  interpolation  in  the  stored  discrete 
distribution  Gx,.  The  azimuthal  deflection  A0X,  is  selected 
from  a uniform  distribution  (0,2v)  about  the  present  direc- 
tion of  travel.  The  new  direction  of  travel  (0,<p)  at  the 
end  of  the  pathlength  increment  is  found  using  AO^,  A<{>  , 
and  the  old  direction  ( 0 ‘ , <f>  ‘ ) in  an  orthogonal  transformation. 

After  each  step,  the  new  particle  position  z is 
calculated: 

z = z'  + Asn  cos6n  (57) 

where  z‘  is  the  former  position.  The  distance  |z  - zBND I 
to  the  nearest  boundary  is  calculated.  (For  the  emission 
calculation,  this  boundary  is  the  emitting  surface.)  If 
this  direction  is  greater  than  the  particle  residual  range 
R (E^) , then  the  trajectory  is  terminated. 

If  the  particle  has  crossed  a boundary  during  a 
step,  then  the  energy  of  the  particle  at  the  boundary  is 
found  by  linear  interpolation.  If  is  the  energy  at  the 
start  of  the  step,  and  <5s  is  the  distance  along  the  particle 
path  from  the  start  of  the  step  to  the  intersection  of  the 
path  with  the  boundary,  then  the  energy  of  the  particle  at 
the  boundary  is 
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6s 


(58) 


ebnd  en 


As 


AE 


N 


N 


The  charge  current  at  the  boundary  is  incremented  by  the 
particle  weight  W.  The  energy  current  is  incremented  by 
W- E~..~  . If  the  boundary  is  the  emission  surface,  the 
particle  weight  W is  scored  in  each  of  the  appropriate 
energy,  9-angle  and  $-angle  bins  of  the  emission  distri- 
bution . 

If  the  boundary  crossed  by  the  particle  is  an  inter- 
face between  different  materials,  the  particle  is  backed 
up  to  the  point  of  intersection  of  the  path  with  the 
boundary.  The  particle  then  starts  into  the  new  material 
with  energy  EBND  and  the  same  direction  of  travel  as  when 
it  reached  the  boundary.  The  first  step  in  the  new  material 
is  a partial  step,  following  the  same  partial  step  procedure 
as  at  the  start  of  the  trajectory. 

The  particle  trajectory  is  terminated  with  the  re- 
sidual range  is  less  than  the  distance  to  the  nearest  bound- 
ary, or  when  che  energy  has  decreased  to  1 keV.  A new 
trajectory  is  then  started.  This  process  is  continued  until 
the  desired  number  of  particle  histories  has  been  calculated. 
In  order  to  preclude  the  disaster  of  exceeding  the  CP  time 
limit  before  obtaining  results,  POEM  checks  on  CP  time  during 
the  random  walk  calculation.  If  the  time  limit  is  approached, 
POEM  will  terminate  the  random  walk  calculation  and  move  on 
to  the  calculation  of  the  results. 

3.3  TRANSITION  ZONE  DOSE  AND  CURRENT  PROFILES 

Extension  of  the  above  electron  source  and  transport 
computational  procedures  to  the  calculation  of  the  dose  and 
current  profiles  near  a material  interface  could  be  straight- 
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forward:  calculate  the  source  densities  and  distributions 

on  both  sides  of  the  interface,  and  calculate  the  transport 
in  the  neighborhood  of  the  interface.  We,  however,  have 
found  it  advantageous  to  be  somewhat  devious.  In  the 
straightforward  approach  we  would  spend  considerable  time 
tracking  particles  in  the  transition  zone  whose  trajectories 
never  cross  the  interface.  The  trajectories  of  these  parti- 
cles are  independent  of  the  nature  of  the  material  on  the 
opposite  side  of  the  interface.  In  particular,  their  tra- 
jectories would  be  the  same  if  the  transport  region  were 
in  electron  equilibrium  (i.e.,  the  region  were  in  an  un- 
bounded medium) . We  know  the  equilibrium  dose  and  current 
(the  kerma  and  the  PC  current  in  an  unbounded  medium) . 

These  are  readily  obtainable  from  tabulations  or  by  analy- 
19  20 

tical  codes.  When  available,  analytical  methods  are 

always  preferable  to  Monte  Carlo  methods.  By  restricting 
the  Monte  Carlo  calculation  to  the  nonequilibrium  contri- 
butions to  the  dose  and  current,  considerable  reduction  in 
variance  is  achieved.  Furthermore,  it  increases  physical 
insight  into  processes  determining  the  characteristics  of 
the  dose  and  current  profiles. 


Consider  a planar  interface  between  two  materials, 

A and  B.  We  wish  to  calculate  the  photo-Compton  (PC)  cur- 
rent and  the  dose  distribution  in  B near  the  interface. 

Let  J be  the  PC  current  in  B of  electrons  arising  in  A, 
and  let  be  the  PC  current  of  electrons  arising  in  B 

which  have  been  backscattered  from  A.  Now  consider  the 
interface  when  the  slab  A is  replaced  by  a slab  B'  of 
composition  identical  to  material  B.  (The  interface  is 
now  a plane  in  the  midst  of  an  unbounded  region  of  material 
B. ) Let  JB,  be  the  PC  current  in  B of  electrons  arising 
in  B'  and  be  the  PC  current  in  B of  electrons  arising 

in  B which  have  been  backscattered  from  B'.  Finally,  let 
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be  the  equilibrium  PC  current  in  B (the  PC  yield  in  an 
unbounded  region  of  B) . Then  the  total  PC  current  in  B 
near  the  interface  with  A is  given  by 


i 


i 

! 


i 


J(z)  = JEQ  + [JA(z)  - JB,(z)] 


+ [JB/A(Z)  " JB/B'(Z)1* 


(59) 


The  first  bracketed  term  in  Equation  (59)  represents 
the  nonequilibrium  component  of  the  PC  current  due  to  emis- 
sion from  the  interfacing  material.  The  second  bracketed 
term  represents  the  nonequilibrium  component  due  to  back- 
scatter  from  the  interfacing  material.  A similar  expression 
can  be  obtained  for  the  dose: 


D(z)  = Deq  + (Da(z)  - (Db, (z) ] 
+ Idb/a(z)  - db/b,(z)] 


(60) 


where  D__-  is  the  kerma. 

By 

Monte  Carlo  calculations  are  performed  to  obtain 
each  of  the  nonequilibrium  terms  for  the  PC  current  - 

JA'  jb'  ' JB/A'  an<^  JB/B'  E<3uat-i°n  (59)  " and  the  cor- 
responding PC  vector  energy  currents.  The  nonequilibrium 

terms  for  the  dose  - DA,  Dg,,  an<^  DB/B'  ^ n EcJuation 

(60)  - are  obtained  by  differentiating  the  electron  energy 

currents.  The  equilibrium  PC  current  can  be  obtained  from 

20 

the  Dellm  and  MacCallum  handbook.  The  kerma  can  be  ob- 

19 

tained  from  standard  tables. 

The  relative  magnitudes  of  the  terms  in  Equations 
(59)  and  (60)  vary  with  photon  energy.  For  example,  for 
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dose  in  a low-Z  material  near  a hiyh-Z  material,  D^fx)  is 
very  large  compared  with  any  of  the  other  terms  for  photon 
energies  less  than  about  100  keV.  But  for  photons  with 
energies  greater  than  about  1000  keV,  incident  through  the 
low-Z  material,  the  t,erm  ( x)  dominates.  For  maximum 

efficiency,  the  number  of  particle  histories  allocated  to 
the  calculation  of  each  of  the  terms  should  be  directly 
proportional  to  its  relative  contribution  to  the  sum. 

It  is  planned  to  automate  this  allocation  procedure.  In 
the  meanwhile,  some  user  intuition  is  required,  perhaps 
aided  by  the  forthcoming  report  on  transition  zone  dose 
and  current. 

POEM  also  calculates  analytic  functional  fits  to 
the  Monte  Carlo  calculated  dose  and  current  profiles.  The 
interest  in  fits  is  twofold.  First,  analytic  expressions 
are  of  course  more  convenient  for  representing  results  than 
are  tabulations  of  Monte  Carlo  data.  Second,  it  is  hoped 
to  "damp  out1'  statistical  fluctuations  in  the  charge  and 
energy  current  profiles.  The  Monte  Carlo  predictions  of 
the  charge  deposition  and  dose  are  obtained  through  first- 
order  finite  differencing  of  the  charge  and  energy  currents. 
For  deposition  zone  number  n. 


Qn  = 


j(z  ■>)  “ j(z  ) 


D = 
n 


- f(zn) 


where  j(zn)  and  fCz^)  are  t^e  charge  and  energy  currents, 
respectively,  at  zone  boundary  z^ . This  is  the  crudest  form 
of  numerical  differentiation.  Even  if  the  statistical  errors 
for  j and  f are  small,  the  errors  foj.  Q and  D will  be  large 
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if  A z is  small.  If  one  could  accurately  represent  j (z)  and 
f(z)  with  analytic  expressions,  then  the  differentiation 
could  be  performed  analytically. 

POEM  fits  the  Monte  Carlo  calculated  charge  and 
energy  current  profiles  to  expressions  of  the  form 

f ( z ; = A exp(Bz  + Cz2  + Dz3).  (63) 

Fits  to  the  deposition  profiles  are  obtained  through  analy- 
tical differentiation  of  the  expression  in  Equation  (63): 

= (b  + 2Cz  + 3Dz 2 ) f (z)  (64) 

dz 

Representative  fits  to  the  charge  current  profile 
and  to  the  energy  deposition  profile  are  shown  in  Figure  11. 
The  fits  to  the  current  profiles  are  consistently  quite  good. 
The  fits  to  the  deposition  profiles,  however,  are  not  uni- 
formly satisfying,  for  example  as  shown  in  Figure  12.  The 
reason  for  the  discrepancy  is  that  for  the  contribution  due 
to  electron  emission  from  the  high-Z  material  - the  D^(z) 
term  - the  various  source  electron  groups  - Compton,  K-photo, 
L-photo,  and  Auger  - produce  characteristically  different 
energy  deposition  profiles  in  the  low-Z  material.  The  fit 
can  be  considerably  improved  by  separately  calculating  and 
fitting  the  profiles  for  each  of  the  source  groups,  and  then 
summing  the  contributions  (Figure  13) . The  POEM  code  has 
this  capability,  which  may  be  exercised  as  an  option. 
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'ZT' 


Figure  11 

RELATIVE  DOSE  PROFILE 


(Plots  courtesy  of  J.C.  Garth,  Air  Force  Cambridge  Research 
Laboratories) 
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ENERGY  DEPOSITION  PROFILE  IN  SILICON 


Section  4 


CROSS  SECTIONS 


ELECTRON  PRODUCTION 


Three  electron  production  mechanisms  are  included 
in  POEM:  the  Compton,  photoelectric,  and  Auger  effects. 

Pair  production  and  secondary  electron  production  are  not 
included . 


4.1.1  Compton  Electrons 

The  density  and  energy-angular  probability  distri- 
bution of  the  Compton  electrons  are  calculated  using  the 
theory  of  Klein  and  Nishina.  The  mathematical  expression 

for  the  total  interaction  cross  section  is  taken  from 

_ 21 
Evans : 


2irr 2 
o 


— ^ 2,(W  - - In  (l+2a) 

cc2  l+2a  a 


+ 2a  In  (l+2a)  - (14.2ci)?  j 


where  a = hw/mc2  and  r = e^/mc2.  The  expression  for  the 

° 30 

energy  distribution  is  taken  from  Evans: 


dp (E ) _ 
dE 


I 2 +(  E V. 

1 hv-E 

2 nv-El  ) 

j \hv-E/ 

a2  hv 

a E | 

where  A is  a constant. 
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l , 

[ Once  the  electron  energy  E is  determined,  the  scat- 

["  tering  angle  0 of  the  electron  is  found  from  the  relativistic 

| expressions  for  conservation  of  energy  and  momentum: 


tan0 


1 

1+ot 


- 1)  - 1 


(67) 


i 


> 

f 


4.1.2  Photoelectrons 

The  total  photoelectric  cross  section  is  found  using 

2 3 

log-log  interpolation  with  the  ENDF/B  photon  cross  section 
compilation.  We  have  extracted  the  photoelectric  cross  sec- 
tions from  the  ENDF/B  compilation  for  the  elements  Z = 1 
through  83,  86,  90,  92,  and  94  and  stored  the  cross  sections 
on  the  photon  interaction  data  tape.  The  tape  also  contains 
the  relative  shell  ionization  probabilities  and  the  K,  mean 
L and  mean  M absorbtion  edges  as  discussed  in  Section  3. 

*"■  e nean  N edge  is  approximated  by  zero. 

The  angular  distribution  of  the  photoelectrons  is 

22 

jalculated  using  the  Fischer  formula  for  low  energy  elect- 
rons, E<Eo: 


dp ( 0 ) _ , sin  30 
de  (1-6 ' cos9) 4 


(68) 


3’ 


3 

l+hv/2mc ?- 


(69) 


and  the  Sauter  formula  for  high  energy  electrons,  E>E  : 

G 


dp  (6)  _ . , si.n36(l  * g"cos0) 
d9  (1  - 3cos0) 4 


3"  = 2 — , T = E/m c2 

T(l+T-  ) ” 1 


(71) 


8.1 


where  the  transition  electron  energy  E is  that  recommended 

24  ° 

by  Dellin  and  MacCallum: 

E = 2.5Z  - 0 . 5E„  (72) 

O K 


4.1.3  Auger  Electrons 

The  K and  mean  L fluorescence  yields  have  been  ob- 

25 

tained  from  the  compilation  of  Bambynek , et  al.  and  from 

19 

the  fitted  data  of  Storm  and  Israel,  where  not  available 
from  the  Bambynek  compilation,  and  stored  on  the  photon 
interaction  data  tape.  The  mean  M fluorescence  yield  is 
approximated  by  zero.  The  angular  distribution  of  the 
Auger  electrons  is  isotropic. 


4.2  ELECTRON  TRANSPORT 

4.2.1  Energy  Loss 

POEM  uses  the  continuous  slowing  down  approximation. 
Energy  loss  AE  is  uniquely  determined  by  pathlength  incre- 
ment As.  The  energy  loss  -dE/ds  includes  two  terms,  colli- 
sion and  radiation: 


dE 

ds 


(73) 


The  collision  term  is  calculated  usi..g  the  expres- 

2 6 

sions  of  Berger  and  Seltzer  for  the  Rohrlich  and  Carlson 

formulation  of  the  Bethe  stopping  power  theory  with  the 

Sternheimer  density  effect  correction.  For  Z>_13,  the  em- 

2 c 

pirical  relation  of  Sternheimer*  is  used  for  the  mean  ex- 
citation energy.  For  Z<13,  the  tabulated  values  of 
Reference  26  are  used. 
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For  the  POEM  code  range  of  application,  the  radia- 
tion loss  term  is  small  compared  with  the  collision  loss 
cerm.  A simple,  approximate  formulation  of  the  radiation 
loss  term  has  proven  adequate.  The  collision  loss  is  cal- 
culated as 


(-S), 


= NE  <i>  (Z)R(E,Z) 


where  N is  the  density  (atoms/cm),  E is  the  electron  energy, 
<i>  is  given  by 


* = 2U+11  r2 

137  o 


where  r = e2/mc2.  R is  given  by  a least-squares  fit  to 
O 2 6 

the  data  of  Berger  and  Seltzer:  for  E<1  MeV,  R is  con- 

stant; for  1 MeV<E<10  MeV,  log  R is  a linear  function  of 
log  E . 


4.2.2 


Scatterim 


The  probability  distribution  GN(0)  for  the  angle  0 

through  which  the  electron  is  scattered  in  traversing  a 

pathlength  increment  AsN  is  calculated  using  the  multiple 

scattering  theory  of  Goudsmitt  and  Saunderson  as  implemented 

3 

by  Berger  in  ETRAN.  The  cumulative  probability  distribution 
is  expanded  in  Legendre  polynomials: 


G.,  ( 6 ) =11  + 

N max  , 


(,t  + 2>  <H*  - "max’ 


P()  (cosO  ' ) sinO  ' dO  1 
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where 


H a = 


I 7“  I 

exp  { - I S£ (s 1 ) ds ' ( 


I 


s 


(77) 


and 


/ 


S£(s)  = 2rN  I o(6,s)  (1  - P?  (cos0) ] d (cos0)  (78) 

-1 


a is  the  Mott  single  scattering  cross  section  with  the 
screening  correction  of  Moliere.  The  coefficients  S p are 
calculated  using  the  procedures  of  the  SING  subroutine  of 
ETRAN.  The  coefficients  H£  are  calculated  using  the  pro- 
cedures of  the  MULT  subroutine  of  ETRAN.  The  term  H 

max 

represents  the  probability  of  zero  scatter  over  the 
pathlength  increment, 

H = Lim  . (79) 

max  £_^rj0  % 


Separation  of  this  term  from  the  Legendre  expansion  consid- 
erably improves  the  series  convergence  at  low  energy,  as 
discussed  in  Reference  12.  The  integrals  of  the  Legendre 
polynomials , 


(cosQ ' ) 


sinO ' 


dO' 


are  evaluated  through  recursion. 


(80) 
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Section  5 


USER  INFORMATION 


POEM  is  a FORTRAN  computer  program  coded  for  the 
CDC  6400,  6600,  or  7600  computer.  (A  version  coded  for 
the  IBM  370/195  computer  is  also  available.)  The  code  is 
obtainable  on  seven  track,  556  bpi  tapes.  Two  tapes  are 
required  for  operation:  one  contains  the  photon  inter- 

action data  file;  the  other  contains  the  program  source 
decks  m the  CDC  UPDATE  format.  Five  codes  are  contained 
on  the  program  tape,  identified  by  UPDATE  DECKNAMES: 

EDITOR,  DATAGEN,  1S1D,  1S2D,  and  2S.  EDITOR  reads  the 
photon  interaction  data  tape  and  writes  a selected  subset 
of  data  on  disk  or  tape.  DATAGEN  creates  an  electron  in- 
teraction data  file  on  disk  or  tape.  1S1D  and  1S2D  are  the 
one-  and  two-dimensional  x-ray  photoemission  transport  pro- 
grams. 2S  is  the  x-ray  transition  zone  dose  transport 
program. 

A transport  calculation  requires  a photon  inter- 
action data  file  (lfn*  = TAPE 8 ) , an  electron  interaction 
data  file  (lfn  = TAPE16)  and  a card  input  which  defines  the 
transport  media,  the  photon  spectrum,  and  electron  source 
parameters.  The  procedures  for  creating  the  data  files 
and  running  the  transport  programs  are  outlined  in  the 
succeeding  paragraphs. 


* local  file  name 


5.1  PHOTON  INTERACTION  DATA  FILE  CREATION 

The  EDITOR  code  permits  the  creation  of  a truncated 
photon  data  set  on  permanent  disk  file  from  a subset  of  the 
photon  data  set  on  the  tape  file.  This  subset  can  be  larger 
than  the  set  required  for  a given  transport  calculation;  the 
transport  program  searches  the  file  for  the  data  required 
for  a particular  _un.  The  recommended  procedure  is  to  create 
a permanent  disk  file  with  as  many  elements  as  anticipated 
to  be  required  for  a series  of  calculations. 

A run  requires  the  EDITOR  program,  the  photon  inter- 
action data  tape  (lfn  = TAPE1) , and  a card  input.  The  pro- 
gram creates  a photon  interaction  data  file  (lfn  = TAPE2 ) 
to  be  cataloged  on  permanent  disk  file. 

The  card  input  consists  of  two  to  four  cards: 


CARD 

FORMAT 

DESCRIPTION 

1 

15 

Number  of  elements  to  be  put 
in  file  (max  25) 

2-4 

1015 

List  of  element  atomic  numbers 
in  ascending  order 

It  is  essential  that  the  element  atomic  numbers  be  listed 
in  ascending  order  so  that  the  transport  program  can  prop- 
erly search  through  the  file. 

The  printer  output  confirms  the  copying  of  the  in- 
teraction data  for  each  element. 

5.2  ELECTRON  INTERACTION  DATA  FILE  CREATION 

The  DATAGEN  code  creates  the  electron  interaction 
data  file  on  disk  or  tape.  A run  requires  the  DATAGEN 
program  and  a card  input.  The  program  creates  a data  file 
(lfn  = TAPE16) . Once  an  electron  data  set  for  a material 

86 


has  been  created  on  permanent  disk  file  or  tape,  that  set 
can  be  used  for  any  transport  calculation  in  that  material, 
for  any  electrons  of  energy  equal  to  or  less  than  the  maxi- 
mum electron  energy  in  the  file  (EEMAX) . 

5.2.1  Input 

The  card  input  is  in  the  NAMELIST  format.  This  for- 
mat is  convenient  for  the  input  of  optional  parameters.  One 
or  two  NAMELISTS  are  used:  POET  and  ELISTl.  EL1ST1  is  re- 

quired only  if  it  is  desired  to  vary  the  logarithmic  spacing 
of  the  transport  energy  grid  as  discussed  in  Section  3.2.3. 
(POET  is  always  required  to  create  a POEM.)  The  NAMELIST 
input  variables  are  defined  below. 

NAMELIST/POET/ 

For  a material  composed  of  a single  element  the 
variables  defining  the  material  are 

Z(l)*=  atomic  number 
A ( 1 ) = atomic  weight 

For  a material  composed  of  several  elements,  the  variables 
defining  the  material  are: 

LMENT  = number  of  elements  (max  8) 

Z ( 1 ) = list  of  atomic  numbers 
A ( 1 ) = list  of  atomic  weights 
W(l)  = list  of  weight  fractions 


* Notation:  For  the  NAMELIST  input  of  an  array,  the  CDC 

computer  demands  the  index  of  the  first  element,  for  example 
Z ( 1 ) =1,  6,  7,  8.  We  denote  herein  all  arrays  by  the  for- 
mat VAR ( 1 ) = . . . . 
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The  variables  defining  the  transport  energy  grid 

are: 

EEMAX  = maximum  electron  energy  in  MeV,  (£10  MeV) 

(If  one  is  creating  an  electron  data  set  on  permanent  disk 
file  or  tape,  then  one  should  set  EEMAX  = 1 0 so  that  the 
file  will  be  valid  for  any  POEM  calculation  for  that  material, 
unless  of  course  one  does  not  anticipate  any  calculation  for 
electron  energies  greater  than  some  limit  EEMAX<10  MeV.) 

NSUB  = number  of  substeps  per  logarithmic  step 
(Equation  (45)] 

If  one  does  not  want  to  vary  the  logarithmic  spacing  of  the 
transport  energy  grid  between  the  maximum  energy  (EEMAX) 
and  the  minimum  energy  (0.001  MeV),  then  one  inputs 

NCYC  = the  number  of  logarithmic  steps  to  reduce 

the  electron  energy  by  a factor  of  one-half 
[Equation  (44)].  The  recommended  value  is  8. 

If  one  does  want  to  vary  the  logarithmic  spacing  of  the 
transport  energy  grid,  as  discussed  in  Section  3.2.3,  then 
one  does  not  input  NCYC  but  inputs  the  NAMELIST  ELISTl  as 
described  below. 

NAMELIST/EL IST1/ 

ELISTl  is  input  only  if  one  wants  to  vary  the  log- 
arithmic spacing  of  the  transport  energy  grid.  The  vari- 
ables are 

NCMAX  = number  of  energy  grid  points  at  which  the 
logarithmic  spacing  is  to  change  (-"ax  10) 

TC(1)  = list  of  energies  (MeV)  at  which  the  spac- 
ing is  to  shift 
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(Unless  one  wants  to  terminate  the  electron  trajectories  at 
some  minimum  energy  greater  than  0.001  MeV,  set  TC(NCMAX)  = 
0.001.) 

NCYCl(l)  = list  of  NCYC  values  (Equation  (44)] 

(NCYCl(l)  determines  the  logarithmic  spacing  between  EEMAX 
and  TC(1);  NCYC1(2)  determines  the  logarithmic  spacing  be- 
tween TC ( 1)  and  TC(2);  ...;  NCYCl(NCMAX)  determines  the 
logarithmic  spacing  between  TC(NMAX-l)  and  TC(NMAX).  See 
Section  3.2.3  for  discussion  of  the  spacing.) 

NDMAX  = number  of  energy  grid  points  at  which  the 
multiple  scattering  distribution  is  to  be 
printed  (max  10) 

ND ( 1 ) = list  of  energy  grid  point  index  numbers  at 
which  the  multiple  scattering  distribution 
to  be  printed  (e.g.,  if  ND(1)  = 1,  10,  25, 
then  the  distributions  are  printed  at  E, 

E i o / and  E2  5 

(It  is  recommended  that  the  multiple  scattering  distribution 
be  printed  at  least  at  every  grid  point  at  which  the  loga- 
rithmic spacing  changes.  If  one  uses  the  option  of  input- 
ting NCYC  in  the  POET  NAMELIST,  then  the  ND  values  are  in- 
ternally set  to  ND ( 1 ) = 1,  NCYC+1 , 2NCYC+1,  ) 

5.2.2  Output 

The  NAMELIST  inputs  are  printed  out,  complete  with 
the  default  values  for  the  variables  not  input. 

The  multiple  scattering  distribution  is  printed  at 
each  of  the  grid  points  ND(1),  ND(2),  . ..,  NDMAX.  This  out- 
put includes: 
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N = energy  grid  point 

E = energy  within  the  energy  bin  for 
which  the  printed  distribution  is 
calculated 

S = residual  range  (g/cm2) 

DELTAS  = pathlength  increment  (g/cm2  and 
in  mean-free-pathlengths) 

MFP  = transport  mean-f ree-pathlength 
(g/cm2 ) 

MEAN (COS)  = mean  cosine  of  angle  through 
which  electron  is  deflected 

PROB (NO  SCATTER)  = probability  of  zero  scatter  along 

pathlength  increment 

SINGLE  SCATTERING  LEGENDRE  COEFFICIENTS 

= r S,, , 9,= 1,  l , where  r =resid- 
o t max  o 

ual  range  and  the  coefficients 

S^  are  given  by  [Equation  (78)] 

MULTIPLE  SCATTERING  LEGENDRE  COEFFICIENTS 

‘ V!W'  l=l'  whare  the 

coefficients  H&  are  given  by 
[Equation  ( 77) ] 


(The  multiple  scattering  coefficients  must  monotonically 
converge  to  zero.) 

CUMULATIVE  MULTIPLE  SCATTERING  DISTRIBUTION 

= G^  ( 0 ) , 0=0,  it,  [Equation  (76)] 

(The  distribution  must  converge  monotonically  to  one.) 

The  end  point  deflection  angles  for  the  multiple 
scattering  distribution  are  printed. 

The  electron  slowing  down  spectrum  is  printed. 

This  output  includes:  N = index  number  of  logarithmic 

grid,  energy,  stopping  power  - collision,  radiation,  and 
total  - step  size  (pathlength  increment) , and  mean  penetra- 
tion in  the  electron  initial  direction. 
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5.3  X-RAY  PHOTOEMISSION  CALCULATION 

The  Monte  Carlo  codes,  deck  names  1S1D  and  1S2D, 
calculate  the  x-ray  photoemission  from  a planar  surface. 

For  cases  of  normal  photon  incidence  or  when  the  two- 
dimensional  description  of  the  emission  angular  distribu- 
tion is  not  required  the  ID  version  should  be  used  since 
it  is  a significantly  faster  code.  A run  requires  the 
transport  program,  a photon  interaction  data  file  (lfn  = 
TAPE8)  containing  the  elements  of  the  material,  the  electron 
interaction  data  file  (lfn  = TAPE16)  for  the  material,  and 
a card  input. 

5.3.1  Input 

The  card  input  consists  of  three  NAMELISTS  - MEDIA 
defining  the  material,  PHOTON  defining  the  photon  spectrum, 
and  ESOURC  defining  source  parameters  - and  a time  limit 
card.  The  deck  structure  is 

MEDIA 

PHOTON 

ESOURC 

time  limit  (CP  seconds  in  F10.0  format) 

The  variables  of  the  NAMELISTS  are  defined  in  the  succeed- 
ing paragraphs. 

NAMELIST/MEDIA/ 

MEDIA  defines  the  electron  emitting  material.  For 
a material  composed  of  a single  element  the  composition  is 
defined  by  a single  variable 

NZEL(l)  = atomic  number 
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For  a composition  of  several  elements  the  variables  are 

NELMTS  = number  of  elements  (max  10) 

NZEL(l)  = list  of  atomic  numbers 
PCOMP(l)  = list  of  weight  fractions 


For  a thick  slab  (Az^R^)  the  emission  is  independent  of 
thickness;  no  input  requiring  thickness  or  density  is  re- 
quired. For  a thin  foil  (Az<Rg)  the  thickness  may  be  input 
in  one  of  three  ways: 

(a)  TKCMS  = thickness  in  g/cm2 

(b)  TKCMS  = thickness  in  cm 

RHO  = density  (g/cm3) 

(c)  TKINS  = thickness  in  inches 

RHO  ~ density  (g/cm3) 

NAMELIST/PHOTON/ 

PHOTON  defines  the  photon  energy  spectrum  and  angle 
of  :idence.  Three  options  are  available  for  defining  the 


spectrum: 


(a)  monochromatic  spectrum 
TBL(l)  = photon  energy  (keV) 

(b)  arbitrary  spectrum 

NHNU  = number  of  photon  energy  groups  (max  120) 

TBL(l)  = list  of  mid-points  of  the  energy  bins 
(keV) 

WIDTH (1)  = list  of  widths  of  the  energy  bins  (keV) 

Eli ( 1 ) = list  of  the  fractions  of  the  photon 
energy  fluence  in  each  energy  bin 

(c)  black  body  spectrum 

T = spectrum  temperature  (keV) 


The  default  angle  of  incidence  is  zero  (normal  incidence) . 
For  oblique  incidence  input 
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ANG  = angle  of  incidence  with  respect  to  surface 
normal  (degrees) 

For  back  emission  one  can  either  input  ANG>90°,  or  input 
ANG<90°  and  set  NSIGN  - -1  in  the  ESOURC  NAMELIST. 

NAMELIST/ESOURC/ 

ESOURC  sets  parameters  defining  the  electron  source. 
The  only  required  input  is 

IMAX  = number  of  particle  histories  to  be  run 
(the  default  value  is  1000) 

Normally  the  program  prints  the  marginal  source  distribu- 
tions S'(E)  and  S'(0)  but  not  the  joint  distribution  S'(E,9). 
To  suppress  the  print  of  the  marginal  distributions  set 

MARGIN  = 0 

To  print  the  joint  distribution  set 
JOINT  = 1 

One  can  specify  back  emission  by  setting  the  photon 
angle  of  incidence  ANG  with  respect  to  the  surface  outward 
normal  greater  than  90°.  Alternatively  one  can  use  the 
angle  with  respect  to  the  inward  surface  normal  (<90°)  and 
set 

NSIGN  = -1 

For  the  2D  code  the  number  of  azimuthal  (<p)  output 
angle  bins  is  8 by  default.  To  change  the  number  of  <j>  bins 
set 

NPHI  = number  of  $ angle  bins  (max  10) 
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The  maximum  electron  energy  in  the  source  distribu- 
tion is  by  default  equal  to  the  maximum  photon  energy, 
TBL(NHNU).  To  reduce  this  maximum  energy  set 

EEMAX  = maximum  electron  energy  (keV) 

The  number  of  election  energy  groups  and  the  number 
of  electron  angle  groups  m the  source  distribution  are  by 
default  equal  to  20,  the  maximum.  To  reduce  the  number  of 
source  groups  (correspondingly  reducing  the  number  of  groups 
In  the  output  distribution)  set 

NEE  = number  of  energy  groups 

NTH  = number  of  angle  groups 

5.3.2  Output 

The  program  prints  out  the  NAMELIST  inputs,  complete 
with  the  default  values  of  the  variables  not  inputted.  The 
photon  spectrum  representation  is  printed.  Confirmation 
is  printed  of  the  elements  read  from  the  photon  interaction 
data  file. 

The  source  material  definition  is  printed:  compo- 

sition, thickness,  and  density.  The  electron  source  density 
S (electrons/cm2 ) and  the  electron  source  density  distribution 
with  respect  to  photon  energy  S ' (hv) (electrons/cm2 ) are 
printed.  For  MARGIN  = 1 (default  value  = 1)  the  one- 
dimensional, marginal,  electron  source  energy  and  angular 
distributions  S'(E)  and  S'(Q)  are  printed.  For  JOINT  = 1 
(default  value  = 0)  the  two-dimensional,  joint  distribution 
S' (E,6)  is  printed. 

The  program  prints  the  case  - fore  or  back  emission  - 
the  photon  angle  of  incidence,  and  the  mean  photon  energy. 
Transport  diagnostic  parameters  are  printed:  the  number 

of  particle  histories  run  (if  the  run  time-limits,  this 
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will  be  less  than  IMAX) , the  number  of  particle  emissions, 
the  fractional  number  of  emissions  per  history,  the  mean 
probability  of  electron  emission  from  the  source  region 
(since  the  particles  have  various  weights,  this  probability 
differs  from  the  fractional  number  of  emissions  per  history) , 
tb  i mean  number  of  trajectory  steps  per  particle  emission, 
and  the  mean  angular  deflection  of  the  emitted  electrons 
(with  respect  to  their  original  direction) . 

The  electron  yield  is  calculated  in  units  of  elect- 
rons/photon, coul/cal,  and  cal (electron) /cal (photon) . 

The  following  moments  of  the  distribution  of  the 
are  calculated: 

= mean  energy 

-+■ 

= <rd'z>  = mean  direction  cosine  with  respect 
to  z (the  surface  outward  normal) 

~y 

= <Q*x>  = mean  direction  cosine  with  respect 
to  x (the  projection  of  the  photon  direction 
vector  on  the  material  surface) 


= yield  (electrons/photon)  = the  z-component 
of  the  emission  current  (normal  to  the 
surface) 

- the  x-componcnt  of  the  emission  current 
(parallel  to  the  surface) 


electron  yield 

<Es 
<cos b> 

<cos20> 

<sinOcosb> 


<cos4-> 

J 


J /J 
x'  z 


The  one-dimensional  distributions  of  the  electron 
yield  are  printed: 

y1 (E)  = energy  spectrum 

y' (0)  = polar  angular  distribution  (with  respect 
to  surface  normal) 
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y' ($)  = azimuthal  angular  distribution  (with  respect 
to  the  plane  of  the  photon  direction  vector 
and  the  surface  normal) (for  oblique  incidence 
only;  for  normal  incidence  y'(4>)  is  uniform 
over  ( 0 , 2tt  ) 


[In  future  versions  of  the  code  it  is  planned  to  add  the  one- 
dimensional distributions 

y ' (cosO ) 
y' (sin6cos$)  ] 

The  two-  and  three-dimensional  distributions  of  the 
electron  yield  are  written  on  disc  for  economy  of  printed 
output  (lfn  = TAPE13) . The  file  may  be  saved  for  later 
analysis  by  cataloging  as  a permanent  disk  file  or  tape  or 
printed  by  rewinding  TAPE13  and  copying  it  into  OUTPUT. 

The  distributions  are 

y'  (e,4>) 
y' (E,e) 
y'  (E,*) 
y'  (E,6,<j>) 


5.4  X-RAY  TRANSITION  ZONE  DOSE  CALCULATION 

Calculation  of  the  transition  zone  dose  and  current 
near  a planar  interface  requires  the  transport  program  - 
deck  name  2S  - a photon  interaction  data  file  (lfn  = TAPE8) 
containing  the  elements  in  the  electron  source  material, 
an  electron  interaction  data  file  for  the  two  interfacing 
materials,  and  card  input.  The  electron  data  file  can  be 
created  in  one  of  two  ways:  (1)  set  MEDMX  = 2 in  the  POET 

NAMELIST  input  for  DATAGEN  for  the  data  calculation  for 
material  A,  and  input  a second  complete  card  input  (without 
specifying  MEDMX)  for  material  D in  the  same  run;  or  (2) 
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independently  create  two  data  files  - one  for  material  A, 
or.e  for  material  B - and  then  create  a single  file  through 
the  procedure: 

COPYBF (A,TAPE16) 

BKSP (TAPE16) 

COPYBF (B,TAPE16) 

REWIND (TAPE1 6) 

5.4.1  Input 

The  same  three  NAMELIST  inputs  and  time  limit  card 
are  required  as  in  the  x-ray  photoemission  calculation  - 
the  deck  structure  is  the  same.  The  PHOTON  and  the  ESOURC 
NAMELIST  inputs  are  the  same.  The  only  differences  are  in 
the  NAMELIST  MEDIA.  Refer  to  Section  3.3  for  the  definition 
of  the  problem  geometry  and  discussion  of  the  computational 
procedures . 

NAMELIST/MEDIA/ 

The  identification  of  the  electron  source  material 
is  controlled  by  MSOURC.  For  MSOURC  = 1 (the  default 
value)  the  source  is  in  material  A;  the  program  calculates 
dose  and  current  distribution  in  material  B due  to  electron 
emission  from  material  A.  For  MSOURC  = 2 the  source  is  in 
material  R;  the  program  calculates  the  dose  and  current  in 
material  B due  to  electron  backscatter  from  material  A. 

The  variables  defining  the  electron  source  material 
composition  are  the  same  as  for  the  x-.-ay  photoemission  cal- 
culation. One  only  inputs  the  composition  definition  vari- 
ables (NELMTS,  NZEL,  PCOMP)  for  the  source  material  - if 
MSOURC  = 1,  material  A;  if  MSOURC  = 2,  material  B. 

If  the  materials  are  each  thicker  than  the  maximum 
electron  range,  then  the  dose  and  current  profiles  are 
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independent  of  thickness.  The  thicknesses  are  by  default 
set  equal  to  the  maximum  electron  ranges  - no  input  is 
required.  For  foils  thinner  than  the  electron  range,  the 
thicknesses  may  be  defined  in  one  of  three  ways: 

(a)  TKCMS(l)  = thickness  (A),  thickness  (B) 

(g/cm2 ) 

(b)  TKCMS(l)  = thickness  (A),  thickness  (B) 

(cm) 

DENS ( 1 ) = density  (A),  density  (B)(g/cm3) 

(c)  TKINS(l)  = thickness  (A),  thickness  (B) 

( inches ) 

DENS ( 1 ) = density  (A),  density  (B)(g/cm3) 

5.4.2  Output 

The  first  part  of  the  program  output  - NAMELIST 
inputs,  photon  spectrum  representation,  elements  read  from 
photon  data  file,  source  material  definition,  and  electron 
source  density  and  distributions  - are  identical  to  the 
x-ray  photoemission  output.  A diagnostic  indicating  the 
number  of  materials  for  which  electron  interaction  files 
were  read  is  printed. 

For  MSOURC  = 1 the  transport  diagnostic  parameters 
printed  are:  the  number  of  particle  histories;  the  number 

of  particles  emitted  from  material  A into  material  B;  the 
number  of  particles  backscattered  from  B;  the  electron 
yield  from  A in  olectrons/photon , coul/cal,  and  cal (electron) / 
cal (photon);  and  the  charge  and  energy  backscatter  yields  from 
B (dimensionless) . For  MSOURC  = 2 the  transport  diagnostic 
parameters  printed  are:  the  number  of  particle  histories; 

the  number  of  particles  emitted  from  B into  A;  the  number 
of  particles  backscattered  from  A,  the  electron  yield  from 
B in  electrons/photon,  coul/cal,  and  cal (electron) /cal (photon) ; 
and  the  charge  and  energy  backscatter  yields  from  A (dimen- 
sionless) . 


98 


For  MSOURC  = 1 the  energy  and  angular  distributions 
of  the  electron  yield  from  A,  y'lE)  and  y ' ( 6 ) , are  printed. 
For  MSOURC  = 2 the  energy  and  angular  distributions  of  the 
electron  backscatter  yield  from  A are  printed. 

The  electron  current  (electrons/photon)  and  energy 
current  (cal (electron) /cal (photon) ) are  printed  at  each 
deposition  zone  boundary.  For  MSOURC  -•  1 this  is  the  cur- 
rent emitted  from  material  A;  for  MSOURC  = 2 it  is  the 
current  backscattered  from  A.  These  currents  are  finite 
differenced  to  obtain  the  charge  deposition,  (electrons/g) / 
(photon/cm2),  and  energy  deposition,  (cal /g) / (cal/cm2 ) = 
cm2/g,  in  each  deposition  zone. 

The  electron  charge  and  energy  currents  are  fit  by 
expressions  of  the  form 

y = A exp(Bz  + Cz2  + Dz3) 

The  fit  coefficients  A,  B,  C,  and  D,  the  fit  representations 
of  the  currents  at  each  deposition  zone  boundary,  and  the 
percentage  difference  between  the  fit  and  Monte  Carlo  repre- 
sentations of  the  currents  at  each  zone  boundary  are  printed. 
The  fit  expression  is  analytically  differentiated  to  obtain 
the  fit  representations  of  the  energy  and  charge  depositions. 
The  fit  representations  of  the  depositions  and  the  percentage 
difference  with  respect  to  the  Monte  Carlo  finite  difference 
representations  are  printed  for  each  deposition  zone. 
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